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This  study  is  intended  to  add  to  the  available  technology  for  the 
modeling  of  electromagnetic  scattering  from  random  rough  surfaces.  Two 
primary  contributions  are  made  in  this  thesis.  First,  a  true  joint 
exponential  probability  density  function  is  introduced  as  a  new  model  for 
rough  surface  representation  and  it’s  electromagnetic  scattering 
properties  are  studied.  Second,  a  digital  technique  is  developed  to 
generate  rough  surfaces  using  only  marginal  density  functions.  The 
electromagnetic  scattering  properties  of  these  surfaces  are  then 
found  using  standard  scattering  prediction  techniques. 
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Abstract 

The  purpose  of  this  investigation  is  to  extend  the  capability  of 
statistical  modeling  for  the  problem  of  electromagnetic  scattering  by 
random  rough  surfaces.  A  majority  of  prior  studies  considered  the  well 
known  Gaussian  joint  probability  density  function  (pdf),  partly  because 
of  its  mathematical  convenience.  About  25  years  ago,  a  second  pdf 
called  exponent ial^type  was  introduced  and  has  since  been  used  by 
several  workers  in  the  field.  In  this  thesis  it  is  shown  that  this 
exponential -type  pdf  fails  to  meet  an  essential  statistical  requirement. 
To  remedy  this  situation,  an  exponential  pdf  that  satisfies  all  the 
statistical  requirements  is  derived.  The  scattering  properties  of 
surfaces  modeled  by  this  new  exponential  pdf  are  then  formulated. 

In  the  second  part  of  this  thesis  a  user-oriented  computer  program 
is  developed  to  create  random  surfaces  whose  statistical  description  is 
at  the  will  of  the  user  requiring  only  the  input  of  a  roonovariate  pdf. 
By  using  this  program,  the  average  radar  cross  section  properties  of 
surfaces  modeled  by  the  Gaussian  and  exponential  pdf’s  are  investigated 
further.  Finally,  applications  of  the  work  completed  in  this  thesis  to 
problems  needing  further  investigation  are  outlined. 


ELECTROMAGNETIC  SCATTERING  BY  AN  EXPONENTIALLY 
DISTRIBUTED  ROUGH  SURFACE  WITH  THE  INTRODUCTION 
OF  A  ROUGH  SURFACE  GENERATION  TECHNIQUE 


I .  Introduction 


Extensive  research  has  gone  into  rough  surface  scattering  theory 
where  the  rough  surface  is  modeled  by  a  joint  probability  density 
function  (jpdf).  In  these  models  the  surface  heights  are  random 
variables.  However,  the  research  has  been  limited  to  the  use  of  only 
two  types  of  joint  densities  as  rough  surface  models.  In  addition,  one 
of  these  two  jpdf’s  has  characteristics  which  make  its  use  as  a  surface 
model  inappropriate.  Therefore,  development  of  scattering  properties  of 
surfaces  with  different  jpdf  representation  as  well  as  a  better  overall 
understanding  of  rough  surface  scattering  is  desirable. 


As  noted  by  Pyati,  for  a  joint  probability  density  function  to 
properly  represent  a  surface,  it  must  not  only  satisfy  the  statistical 
properties  required  of  a  second  order  jpdf  but  it  also  must  satisfy  any 
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additional  conditions  imposed  by  the  physics  of  the  problem  (4:3). 
Clearly,  two  such  physical  conditions  should  be  that: 

(1)  As  the  distance  between  any  two  points  ( r)  on  the 
surface  goes  to  infinity,  the  heights,  C  and  ,  at 
the  two  points  should  be  statistically  independent. 

(2)  .As  the  distance  between  any  two  points  on  the  surface 
goes  to  zero,  the  probability  that  the  two  points  are 
equal  in  height  should  go  to  one. 

The  two  jpdf’s  that  are  commonly  used  in  the  current  literature  to 
represent  random  rough  surfaces  are  the  joint  Gaussian  pdf  (Eq  (1)) 


P2(Wr)  = 


2"h2 [  1  -  P2 ( r)] 


1  ✓  2 
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-  2<1<2  P(r>  ♦  <2 


2  h' 


*  - 
1  -  P2  (  r) 


(1) 

2=  surface  heights  at  points  (xt ,y^  )  and  (x2,y2> 

P(r)  s  the  correlation  function 
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h  s  variance  of  ^ 
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and  the  joint  exponential- type  pdf  {Eq  (2)). 
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The  joint  Gaussian  pdf  satisfies  conditions  (1)  and  (2)  and  is  a  valid 
jpdf  for  representing  rough  surfaces .  The  exponential-type  jpdf  fails 
to  satisfy  condition  1,  some  consequences  of  which  will  be  shown  in  this 
thesis,  and  may  therefore  not  be  valid  for  rough  surface  representation. 
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Pyati  has  derived  the  following  true  exponential  joint  probability 


density  function  that  satisfies  both  of  the  physical  conditions  given 
above  [17:8]: 


P2  <Wr)  = 


where 


4h2  (l-p2} 


exp 


s 

U  I  +  U  1 

•  i  •  i  2  1 

h  [1-p2] 

zpc  i<1  i  u2n 


1/2 


h  (i-p2) 


+  0/n2)  sgn(<1)  sgn(<2) 


{2n  +  l) 


2ivM 


n  =  o 


2PC  KJ  K2I) 


1/2 


h  Cl-P2) 


(3) 


th 


I  3  n  order  modified  bessel  function  of  the  first  kind 
n 


sgn(x)  s  the  sign  of  x  or  zero  if  x  =  zero. 


The  scattering  properties  of  surfaces  modeled  by  the  above  exponential 
jpdf  will  be  studied  here. 

In  addition  to  the  lack  of  density  functions  to  use  as  models, 
there  has  heretofore  never  been  presented  a  method  by  which  the  the 
predicted  scattering  from  the  existing  models  could  be  checked  against 
the  scattering  from  real  surfaces  that  meet  the  statistics.  This  thesis 
wi l 1  address  th i s  problem . 
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2.  Statement  of  the  Problem 

The  problem  is  divided  into  two  primary  parts.  Part  A  will  involve 
the  prediction  of  the  average  RCS  per  unit  area  of  a  surface  represented 
by  Pyati’s  exponential  jpdf  using  somewhat  traditional  rough 
surface  scattering  theory.  Part  B  involves  the  development  of  an 
entirely  new  approach  to  the  prediction  of  rough  surface  scattering  for 
surfaces  represented  by  probability  density  functions. 


and  Limitations 


Hie  following  assumptions  will  be  made  in  this  thesis: 

(1)  The  radii  of  curvature  of  the  surface  irregularities  are  large 
with  respect  to  a  wavelength.  This  restriction  permits  one  to 
use  the  physical  optics  method  to  solve  the  electromagnetic 
integral  equations  introduced  in  the  next  chapter. 

(2)  Multiple  scattering  and  shadowing  will  be  ignored. The  effects 
of  shadowing  are  being  accounted  for  in  another  thesis  by  Reif 
[19]  concurrent  with  this  one. 

(3)  The  surfaces  will  be  assumed  to  be  perfect  conductors. 

(4)  The  correlation  length  £  is  assumed  to  be  much  less  than  the 
either  dimension  of  the  surface.  For  instance,  if  the  surface 
is  chosen  to  be  a  rectangle  with  transverse  dimensions  of  Lx 
and  Ly,  £  “  Lx,  Ly.  This  condition  makes  certain  that  the 
scattering  is  a  result  of  the  ensemble  surface  roughness 
rather  than  from  a  statistically  small  portion  of  the  surface 


that  may  not  be  representative  of  the  surface  as  a  whole. 


(5)  The  Gaussian  correlation  coefficient  model, 

2  2 

P(r)  -  exp  [-r  /£  ],  will  be  used  be  used  to  maintain 
consistency  with  previous  calculations  based  upon 
Eqs  (1)  and  (2),  [20:721]  thus  allowing  direct 
comparisons  of  the  results. 

4.  Approach 

For  Part  A  of  the  analysis  a  general  solution  for  the  average  RCS 
per  unit  area  of  rough  surfaces  represented  by  joint  probability  density 
functions  is  formulated  from  the  Stratton-Chu  integral  equation  for 
electromagnetic  scattering.  The  general  solution  arrived  at  is  that 
which  was  previously  obtained  by  Beckmann  [5].  The  general  solution  is 
then  applied  to  Pyati’s  exponential  jpdf  and  the  results  are  compared  to 
those  for  Gaussian  and  exponential-type  surfaces. 

In  Part  B,  a  technique  is  developed  that  allows  the  creation  of 
correlated  rough  surfaces  from  the  one  dimensional  marginal  pdf  of  the 
surface  heights.  This  technique  is  then  used  to  create  rough  surfaces 
that  obey  Gaussian  and  exponential  statistics.  The  physical  optics 
radar  cross  sections  of  these  surfaces  is  then  calculated  using  a 


standard  RCS  prediction  software  package.  These  results  are  then 
compared  to  the  predictions  in  Part  A. 


The  only  support  required  for  Part  A  is  an  IBM  personal  computer 
system  and  associated  software,  all  of  which  are  the  personal  property 
of  the  author.  The  RCS  prediction  package  required  for  Part  B  is  owned 
by  AFWAL/AAWP-3  and  is  hosted  on  their  VAX  11/785.  A  personal  copy  of 
this  code  along  with  access  to  the  VAX  was  provided  by  AFWAL/AAWP-3. 


In  this  section  two  equations  for  the  average  RCS  per  unit  area  of 
random  rough  surfaces  will  be  derived.  The  derivation  of  each  equation 


involves  assumptions  and  approximation  that  will  be  pointed  out  as 
they  occur. 

I .  Stratton-Chu  and  Helmholtz  Equations 

A  starting  point  for  discussions  of  electromagnetic  scattering 
theory  is  the  set  of  equations  known  as  the  Stratton-Chu  integral 
equations  (Eqs  (4)  and  (5))  [23:466]. 
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The  above  equations  were  derived  by  appLying  the  vector  analog  of 
Green’s  theorem  to  Maxwell’s  field  equations  123:464,465).  The  left 
side  of  the  equations  represent  the  complex  electric  and  magnetic  field 
vectors  at  a  point.  P’  located  within  the  closed  volune  of  the  first 


( 

t 
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integral  of  either  equation.  The  closed  surface  of  the  second  integral 
is  that  surface  which  encloses  the  volume  of  the  first  integral.  Within 
the  context  of  these  equations,  p  is  the  charge  density  within  the 
closed  volume  and  J  represents  any  current  sources  within  the  volume. 

The  operator  (*)  indicates  that  the  complex  conjugate  of  the  quantity  is 
to  be  used.  The  E’s  and  H’s  within  the  second  integral  of  Eqs  (4)  and 
(5)  represent  the  total  electric  and  magnetic  fields  at  the  interior 
surface  of  the  volume.  The  symbol  y>  is  the  free- space  Green’s  function. 
This  scenario  is  depicted  in  Fig.  1. 


Fig.  2.  Rough  Surface  Geometry 


The  Stratton-Chu  equations  can  be  simplified  somewhat  since  the 
environment  surrounding  the  target  can  be  assumed  to  be  source  free.- 
This  may  seem  to  be  a  contradiction  since  there  must  obviously  be  some 
source  of  the  original  electromagnetic  fields  that  are  to  be  scattered  by 
the  target.  This  contradiction  is  resolved  by  realizing  that  when  the  RCS 
of  a  target  is  being  considered,  only  those  fields  scattered  from  the 
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target  in  the  direction  of  the  receiver  are  of  any  concern.  Certainly 
there  may  be  additional  fields  at  the  observation  point  P’ .  But,  if 
they  are  not  traveling  along  a  vector  in  a  direction  from  the  target,  or 
are  not  a  result  of  the  scattering  of  the  transmitted  wave  from  the 
target,  they  may  be  set  to  zero.  With  this  simplification  then,  the 
Stratton-Chu  equations  can  be  written  without  the  first  integral  of  each 
equation. 


E(P*  ) 


=  -  O  jju/J  (  n  x  H  ]  V  +  [  n  x  E  )  x  V  tp  +  [  n  •  E  ]  7  ¥>  j- 


ds  ( 6 ) 


H(P’  ) 


E]¥>  -  (n  x  H  )  x  7  y> 


(7) 


Although  the  scattering  geometry  of  the  rough  surfaces  under  study 
do  not  coincide  with  the  geometry  for  which  the  Stratton-Chu  equations 
were  derived  these  equations  can  still  be  used  if  one  simply  considers 
the  rough  surface  to  be  a  portion  of  a  closed  surface  (Fig.  (3)). 
Electromagnetic  fields  will  be  present  at  the  rough  surface  but  will  be 
set  to  zero  over  the  remainder  of  the  imaginary  closed  surface.  Setting 
the  fields  outside  of  the  rough  surface  equal  to  zero  causes  the 
tangential  fields  over  the  entire  closed  surface  to  be  discontinuous. 
Stratton  shows  [23:468-470]  that  the  discontinuity  can  be  accounted  for 
by  the  addition  of  a  line  integral  around  the  contour  of  the  open 
surface.  With  this  addition  the  Stratton-Chu  equations  are  given  by  Eqs 


(  8 )  and  ( 9  )  . 


These  equations  can  be  transformed  into  yet  another  form  [23:469]: 
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The  first  and  last  integrals  of  each  equation  come  from  Eqs  (8)  and  (9) 
by  applying  Green’s  theorem  relating  surface  integrals  to  line 
integrals.  The  second  integral  is  of  course  the  added  term  for  the  edge 
scattering  effects.  The  first  integrals  of  Eqs  (10)  and  (11) 
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are  known  as  the  Helmholtz  Integral  equations.  The  Helmholtz  equations 
are  exact  for  source  free  closed  surfaces. 


2 .  Approximations 

Exact  solutions  to  the  equations  of  the  previous  section  have  not 
been  found  except  in  very  limited  cases.  In  general,  any  solution  to 
these  electromagnetic  scattering  equations  is  facilitated  by  the 
imposition  of  one  or  more  approximations.  This  subsection  discusses 
several  approximations  that  will  be  made  in  this  thesis. 

2 . 1  Far  Field  Approximation  and  Plane  Waves 

If  one  only  considers  observation  points  (P’s)  that  are  very  far 
from  the  scattering  surface  several  simplifications  can  be  made.  The 
first  simplification  is  that  the  incident  fields  at  the  surface  may  be 
considered  plane  waves.  This  will  be  useful  when  the  physical  optics 
approximation  is  made  and  the  angle  of  incidence  at  each  point  along 
the  surface  needs  to  be  known.  If  the  incident  fields  are  not 
considered  to  be  planar,  determining  this  angle  of  incidence  becomes  a 
significant  problem. 

The  incident  planar  electric  field  can  be  written 
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E  =  -  exp 

4*R 
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X  v  Z 

Because  the  field  quantity  is  only  the  real  part  of  Eq  (14),  the  sign  of 
the  remaining  complex  exponential  term  is  inconsequential.  Since  the 
time  dependence  of  the  fields  is  henceforth  going  to  be  suppressed,  the 


negative  sign  will  be  used.  This  will  make  the  remaining  term  in  the 


exponential  positive.  The  complex  electric  and  magnetic  fields  can  now 
be  written: 


E  =  -  exp(jk'r)  ;  H  =  -  exp(jk*r)  (15) 

4  nR  4  tiR 

.Another  important  result  of  the  far  field  approximation  is  that 
terms  containing  higher  orders  of  1/R  may  be  neglected.  Letting  V  be  a 
unit  vector  directed  from  a  point  on  the  scattering  surface  to  the 
observation  point  P’ 


n  = 


(16) 


the  following  approximations  that  will  prove  useful  may  also  be 
obtained  [2:19]. 
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2.2  The  Physical  Optics  Approximation 

The  Stratton-Chu  and  Helmholtz  equations  are  of  a  family  called 
integral  equations.  That  is,  the  term  on  the  left  of  each  equation  is 
equal  to  the  integral  of  some  function  of  itself.  This  is  easier  seen  by 
expanding  the  expressions  for  the  fields  along  the  interior  surfaces  as 


E  ( x,y,z ) 
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while  noting  that  the  left  side  of  each  equation  is 

or  H  (x,y,z) 

p,  scattered  p, 

In  general,  integral  equations  of  this  type  are  not  easily  solved.  In 
fact,  closed  form  solutions  of  the  scattering  integrals  have  been  found 
for  only  a  few  simple  geometries  such  as  spheres  and  infinite 
cylinders  (The  solutions  for  these  shapes  is  facilitated  by  the  fact  that 
their  geometries  are  specified  by  a  single  coordinate  of  an  orthogonal 
coordinate  systems).  Therefore,  several  approximate  solution  methods  to 
the  scattering  problem  have  been  formulated.  One  such  approximation  is 
called  physical  optics.  The  physical  optics  concept  is  to  approximate 
the  fields  at  the  surface  by  some  function  of  the  incident  fields  only. 
Thus,  the  integral  scattering  equations  become  definite  integrals. 

The  most  common  physical  optics  approximation  is  the  tangent  plane 
approximation  where  the  fields  at  each  point  on  the  surface  are  set 
equal  to  that  field  which  would  be  present  at  the  point  if  the  surface 
were  an  infinite  plane  tangent  to  that  point.  "Then  the  incident  field 
at  every  point  on  the  surface  may  be  broken  up  into  components,  E1  in 


E  (x,y,z) 

scattered 


the  plane  of  incidence,  and  .  normal  to  the  plane  of  incidence.  The 
scattered  field  given  by  each  of  these  incident  field  components  is  then 
given  by"  [20:55] 


—9  _  —  1 

E  =  R  E 


sp«c  — s 


where  k 


and  [5:21] 
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is  a  unit  vector  in  the  specular  direction  of  scattering 
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In  the  above  equations  Y  is  the  intrinsic  admittance  of  the  reflecting 
medium,  ©^  is  the  incident  single  with  respect  to  the  surface  normal,  and 
the  propagating  medium  is  taken  to  be  free  space. 

The  tangent  plane  approximation  is  valid  only  when  the  radius  of 
curvature  at  points  along  the  surface  are  much  greater  than  the 
wavelength  of  the  incident  fields  thus  limiting  the  effects  of  the  phase 
shift  caused  by  non-flat  surfaces.  This  will  not  be  discussed  in  detail 
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since  there  are  many  fine  references  that  cover  physical  optics 
[20:55-59],  [10:119-130]  and  [14:454-458].  The  important  consequence  of 
this  fact  is  that  only  rough  surfaces  with  no  sharp  peaks  or  valleys 
will  be  covered  in  the  pursuing  calculations.  The  Gaussian  correlation 
coefficient  used  in  this  thesis  forces  surfaces  to  meet  this  criteria  as 
long  as  the  wavelength  is  kept  small  compared  to  a  correlation  length. 


2.3  The  Infinite  Conductivity  Approximation 
For  any  surface  for  which  z  is  not  constant  (see  Fig.  2)  the 
normal  to  the  surface  is  a  function  of  x  and  y.  This  will  make  the 
plane  of  incidence  and  angle  of  incidence  a  function  of  x  and  y.  This 
in  turn  complicates  the  determination  of  the  local  reflection 
coefficients.  A  common  practice  which  reduces  this  complication  is  to 
assume  that  the  scatterer  is  a  perfect  conductor  so  that  the  admittances 
in  Eqs  (25)  and  (26)  go  to  infinity  and  thus  R^  =  -1  and  =  +1  .  The 
scattered  electric  and  magnetic  fields  are  then  given  by  the 
relationships 
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This  results  in  the  total  fields  at  the  surface  being  given  by 


n  <  E  _  =  n  *  (  E  +  E  )=0 

—  surf  —  1  —  scat 


n  v  H  „  =  2n  x  H 

—  surf  —  1 


This  simplifies  the  problem  in  that  the  angle  of  incidence  dependence  of 
the  surface  fields  is  removed.  The  dependence  on  the  normal  to  the 
surface  at  a  given  point  remains. 


3.  Solution  Ignoring  Surface  Slopes 

In  this  section  an  equation  that  gives  the  average  RCS  per  unit 
area  for  a  random  rough  surface  will  be  developed  by  starting  with  Eqs 
(8)  and  (9).  Applying  Stokes  theorem  to  the  line  integrals  one  finds 
that  [2:18]: 
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By  substituting  Eqs  (17),  (19)  and  (20)  into  Eqs  (31)  and  (32), 


substituting  these  results  into  Eqs  (8)  and  (9)  one  can  write  (2:21) 
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while  Eq  (9)  effectively  becomes 
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where  Z  and  Y  are  the  intrinsic  impedance  and  admittance  of  the 
propagating  medium.  Using  the  infinite  conductivity  identity  that  the 
tangential  component  of  the  surface  electric  field  is  zero  and  that 
n  '<  H  =  2n  *  H  1  ,  these  equations  become 
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Since  either  the  electric  fields  or  the  magnetic  fields  at  the 
observation  point  may  be  used  to  determine  the  RCS  of  the  scattering 
surface,  Eq  (36),  being  the  simpler  of  the  two,  will  be  used. 


Recalling  that  n  is  the  unit  vector  originating  at  a  point  on  the 
surface  and  in  the  direction  of  the  observation  point  P’ ,  the  far  field 


condition  gives  that  this  vector  will  be  constant  over  the  surface  of 
integration.  Therefore  this  cross  product  will  be  taken  outside  of  the 
integrand  using  the  relationship  that  A  '  B  =  -  B  A.  Also,  using  the 
rectangular  coordinate  system  with  the  height  of  the  surface  above  the 
z=0  plane  be  given  by  the  random  variable  <.(\,y),  the  following 
identities  can  be  written: 
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Substitution  of  Eqs  (37)  and  (38)  into  Eq  (36)  results  in  [2:33^ 
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At  this  point  the  approximation  that  the  slope  terms  are 


very  small  will  be  applied  reducing  Eq  (39)  to 
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From  Eq  (15)  the  incident  magnetic  field  can  be  written 


-i 


where  the  polarization  and  time  dependence  of  the  field  are  contained  in 
HQl  and  r  is  given  byr  =xx+yy+<(x,y)  z.  The  magnetic  field  can 
then  be  written  [2:41] 
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Restricting  the  derivation  to  the  case  of  backscatter  where  k1  is  equal 
to  -kh  and  thus  k1  equals  -n ,  Eq  (42)  becomes 
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Equation  (48)  is  easily  checked  for  validity  by  setting  the 
random  surface  height  variable  equal  to  zero.  With  this  substitution, 
the  surface  slopes  are  exactly  zero  and  so  Eq  (48)  should  give  the 
physical  optics  solution  for  a  flat  plate.  Direct  integration  of  Eq 
(48)  with  t(x,y)=0  and  the  application  of  the  RCS  formula  given  by 

4  HR  2  <  j  H  (  P  ’  )  |>2 

<RCS>  =  -  (49) 
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results  in  the  RCS  being 
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Equation  (50)  is  the  physical  optics  solution  for  the  RCS  of  an 
infinitely  conducting  flat  plate  with  sides  and  [10:121]. 

The  next  step  then  is  to  take  the  average  of  Eq  (48)  with  respect 
to  the  random  variable  < ( x , y ) . 
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The  solution  above  gives  the  expected  value  of  H(P’).  This  value  is 
relatively  meaningless  for  very  rough  surfaces  because  the  magnetic 
field  is  a  complex  quantity.  The  large  variations  in  phase  caused  by  the 
rough  surface  will  cause  the  phase  of  this  term  to  vary  wildly.  Therefor 


nnuriv  w  l  iw  w  ltt 


<H<P’)>  gives  little  information  about  <|H(P’)|>  which  is  the  contributor 
in  incoherent  measurement  systems.  To  clarify  this  a  simple  example 
will  be  given. 

The  expected  value  of  a  complex  number  is  written  as 


<H(P’)>  =  <Re  H(P’)>  +  j<Im  H(P’)> 


Table  I  shows  what  misleading  things  can  happen  in  the  extreme 
case  of  180°  phase  variation 


TABLE  I 

Comparison  of  the  Average  of  a  Complex  Quantity 
with  the  Average  of  its  Magnitude 

Re  H(P’)  I  Im  HIP’)  I  <Re  H(P')>  j  <Im  H(P')>  I  <|H(P')|> 


With  this  in  mind,  its  obvious  that  in  order  to  find  the  average  RCS  of 
the  surface,  Eq  (49)  one  must  solve  for  <|H(P’)j>2. 
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The  expected  magnitude  of  the  magnetic  field  at  P’  ,  <  |  H(P’  I  |  >  ,  cai 
be  written  <H( P’ )H( P’  ) *>.  Multiplying  Eq  (51)  by  it’s  complex 
conjugate  results  in 
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By  combining  the  two  integrals  and  taking  the  average  with  respect  to 
the  random  surface  heights,  the  expected  value  of  the  magnitude  of  the 


scattered  magnetic  field  can  be  written 


At  this  time  a  substitution  will  be  made  for  the  sake  of  shorthand 
and  also  to  match  more  closely  the  notation  of  the  common  literature  on 
this  subject  [5],  [20]and  [16].  The  following  parameters  are  defined 
and  simplified  for  the  case  of  backscatter 
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With  this  notation  Eq  (53)  becomes 
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Here  has  been  left  in  for  the  time  being  although  it  is  equal  to 
zero  for  the  case  of  backscatter.  The  inclusion  of  the  limits  of 
integration  indicate  that  each  surface  integral  is  to  be  integrated  over 
a  rectangle  that  is  Lx  by  Ly.  It  will  be  seen  later  that  the  shape  and 
size  of  the  surface  is  inconsequential  so  long  as  it  is  many  correlation 
lengths  in  the  x  and  y  directions. 
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The  averaged  term  within  the  integrand,  <exp( jv  (<-<’))> ,  is  by 
definition  the  joint  characteristic  function  of  the  pdf  evaluated  for  v 
and  -v  [7:428]  and  will  be  denoted  by  X ( v  ,-v  ;r).  It  is  the  two 

z  z  z 

dimensional  fourier  transform  of  the  joint  pdf  and  is  a  function  of  the 

correlation  coefficient  which  is  itself  a  function  of  the  distance 

between  the  two  points  (x,y)  and  (x’,y’)  given  by  r. 

Barrick  [2:45-48]  shows  the  steps  necessary  to  make  a  total 
transformation  of  Eq  (57)  into  a  polar  coordinate  integral  that  car!  be 
integrated  over  the  characteristic  function.  First,  define  the  distance 
between  points  sis 

2  2 

r  =  (x-x')  r  =  ( y— y ’ )  and  r  =  sqrt(  r  +  r  )  (58) 

x  y  x  y 
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The  polar  coordinate  r  has  already  been  defined  as  sqrtfr^ 


2 

+  r 


V 


SO 


all  that  remains  is  to  define  d  and  the  limits  of  integration.  Since 


,(>  .<«  >.l  ,«t 
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neither  of  the  data  points  (x,y)  or  (x’,y  )  are  sis  yet  fixed,  the  easiest 
way  to  define  b  is  as  Barrick  did  [7:79]: 


x-x’  =  r  cos ( d )  and  y-y ’ =  r  sin(b) 


This  definition  is  especially  facilitating  since  the  left  hand  terms 
appesir  in  the  integrand.  Transforming  the  limits  of  integration  is  not 
so  trivial  a  task  but  is  facilitated  by  the  fact  that  for  very  rough 
surfaces,  the  joint  characteristic  function  goes  to  zero  monotonically 
within  only  a  few  surface  correlation  lengths.  Reserving  the 
definitions  of  "very  rough"  aund  "a  few  correlation  lengths"  until  the 
ainalysis  section  where  these  terms  are  examined  more  closely,  suffice  it 
to  say  that  this  allows  the  finite  surface  to  be  expanded  into  an 
infinite  surface  without  affecting  the  value  of  the  integral. 
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The  first  integration  is  trivial  resulting  in 


<|H(P')|>2  =  "  C°V9)  Lx  Ly  |ff/|2 
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2  71  00 
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expj^j(vxr  cos(0)+v^r  sin(0)) 
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2  2 


(62) 


Before  proceeding  the  difference  between  coherent  and  diffuse 
scattering  must  be  understood.  Coherent  scattering  is  that  portion  of 
the  scattered  field  that  follows  the  same  pattern  as  would  the  scattered 
field  if  the  surface  were  not  rough.  For  the  rectangular  surfaces  under 
consideration  here,  the  coherent  pattern  is  that  given  by  Eq  (50).  As 
the  surface  gets  progressively  rougher,  diffuse  scattering  will  become 
more  dominant.  The  diffuse  scattering  is  the  random  scattering  pattern 
that  tends  to  fill  in  the  nulls  of  the  coherent  pattern  and  reduce  the 
power  in  the  coherent  lobes.  A  very  rough  surface  is  one  which  has 
only  diffusely  scattered  power. 

If  a  surface  is  not  very  rough  and  has  a  coherent  component 
remaining  in  the  scattered  field,  the  expansion  of  the  surface  to 
infinite  dimensions  will  cause  the  RCS  to  also  go  to  infinity.  This 
problem  can  be  avoided  if  one  subtracts  out  the  coherent  part  of  the 

power  within  the  integration.  The  coherent  power  can  be  eliminated  by 

—  — t  —  — t 

subtracting  <  H  ><  H  >  from  the  <H  H  > .  This  is  equivalent  to 
subtracting  the  square  of  the  one  dimensional  characteristic  function 


from  the  joint  (two  dimensional)  characteristic  function.  With  this 
modification  and  integrating  over  d  as  shown  by  Beckmann  [5:181-184] 


<x> 
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where  T  =  t(v  ,-v  ;r)  and  t  =  4  (v  )  is  the  characteristic  function  of 

2  z  Z  1  1  z 

the  one  dimensional  marginal  distribution  which  goes  to  zero  as  the 
surface  gets  rougher  and  v  increases.  This  is  easily  converted  to  a 
radar  cross  section  per  unit  area  by  dividing  by  the  area  of  the  surface 
and  employing  Eq(49). 


<RCS> 
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(64) 


Recall  that  the  above  Eqs  (63)  and  (64)  give  only  the  diffusely 
scattered  is  power  and  are  only  valid  for  surfaces  where  the  surface 


slopes  can  be  ignored. 
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4.  Beckmann  s  Solution 

Beckmann  starts  with  the  Helmholtz  integral  equation  for  the 
electric  field  given  by  Eq  (12)  and  repeated  here  [5:19]: 


1  [  f  -  d  f  d  E  ■) 

£<p’>  =  rn  1  r-  - 

J  L  an  d  n  J 


Obviously,  Beckmann  has  ignored  the  line  integrals  introduced  by  the 
fact  that  the  scatterer  is  not  a  closed  surface  .  But,  as  stated  by 
Barrick  states,  these  terms  are  small  and  proportional  to  sin(@)  [2:14]. 
Therefore,  Beckmann’s  solution  should  not  suffer  from  this  deficiency 
at  near  normal  incidence.  As  the  angle  of  incidence  increases,  the 
shadowing  function,  S (9),  begins  to  dominate  the  solution,  again  making 
the  line  integral  terms  inconsequential  (see  Reif  [19]). 

Continuing  with  the  Beckmann’s  derivation,  the  normal  derivative  of 
any  quantity  is  given  by  the  relationship  [22:62] 
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Using  the  7®  given  by  Eq  (20)  and  the  scalar  quantities  for  the 
electric  fields  Beckmann  gives  the  following  physical  optics 
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Continuing  with  Beckmann  s  derivation,  Eq  (68)  is  substituted  into 


Eq  (70)  for  the  surface  normal  vector  and  Eq  (38)  is  substituted  in  for 


the  differential  surface.  Also,  the  magnitude  of  the  Fresnel  reflection 


coefficient  is  replaced  with  unity  for  the  case  of  perfect  conductivity 


yielding 


p  ^  ‘  J  ^  ^  _  ,  A  ,  A  A  _  ^ 

E(P)  =  - -  (k1)*  (—  <xx  —  K^y  +  z  )exp(2jkl*  r)  dxdy  (71! 


Instead  of  simply  dropping  the  slope  terms  in  Eq  (71)  Beckmann 


removes  them  by  employing  integration  by  parts  legitimately.  Expanding 


the  integral  portion  of  Eq  (80)  it  can  be  written  as 


I  =  (-k  <  -  k  K  +  k  )  exp[(2j)(k  x  +  k  y  +  k  ()]  dxdy  (72) 
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For  the  case  of  backscatter  where  k^  =  0,  Eq  (72)  can  be  written  as  the 


sum  of  two  integrals 


(-kx<x>  exp[(2j)(kxx  +  k^K) ]  dxdy 
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Upon  integrating  the  first  integral  by  parts  one  gets 
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Beckmann  shows  that  the  first  integral  represents  an  edge  term  and  is 
negligible  [5:31].  Summing  the  second  and  third  integrals  one  gets 


I  = 
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Since  k  =  ksin(O)  and  k  =  kcos(9),  the  terms  can  be  summed  to  yield 
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Thus 
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E(p)  =  J-^R  -  sec  (0) 


exp[(2j)(k  x+k  £  >]  dxdy 


(77) 


Equations  (77)  and  Eq  (48)  are  practically  identical  in  form  except  the 
cos(&)  term  in  Eq  (48)  is  replaced  by  a  sec(t?)  term  and  the  magnetic  field 


Ms 


is  replaced  by  the  electric  field  in  Eq  ( it ).  Therefore,  the  same 
procedures  can  be  used  to  transform  Eq  (77)  as  were  used  to  transfoim 
Eq  (48)  into  Eq  (64)  and  the  final  result  for  the  average  RCS  per  unit 

area  is 


<RCS>  = 


2k2sec 2 {0) 


rJo(r 
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This  result  is  the  same  as  that  of  Beckmann  [5:79]  with  the  exception 
that  Beckmann’s  result  is  in  terms  of  average  reflectance.  Beckmann’s 
solution  is  also  true  for  the  more  general  case  of  three  dimensional 
bistatic  scattering  [5:720].  Transforming  Beckmann’s  three  dimensional 
bistatic  scattering  with  finite  conductivity  solution  into  RCS  units  one 


<RCS>  = 
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where  I  ■'?  |  is  a  polarization  scattering  matrix  whose  elements  are  given 

in  reference  [5],  For  the  case  of  backscatter  and  infinite 
2  2 

conductivity, I  3  I  reduces  to  sec  (#). 
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5.  Summary 


In  this  section  two  equations  were  derived  to  give  the  average 
diffuse  RCS  per  unit  area  from  a  random  rough  surface.  The  first 
solution  ,  Eq  (64),  was  found  by  making  the  assumption  that  the  surface 
slopes  were  essentially  zero.  The  second  solution,  Eq  (78),  was 
facilitated  by  the  fact  that  the  contour  integrals  around  the  outer  edge 
of  the  surface  were  set  to  zero.  In  comparing  the  Eqs  (64)  and  (78)  it 
is  clear  that  the  only  difference  between  them  is  that  Eq  (64)  has  a 
cosine  term  where  Eq  (78)  has  a  secant  term.  Thus  at  incident  angles 
very  close  to  zero,  these  two  solutions  will  give  the  same  results.  As 
will  be  seen  in  the  next  section,  the  diffuse  backscatter  from  surfaces 
with  small  surface  slopes  is  concentrated  within  small  angles  of 
incidence  and  thus  Eqs  (64)  and  (78)  will  give  essentially  the  same 
results.  .As  the  slopes  of  the  surfaces  increase,  the  diffuse  power  is 
backscatter  is  present  over  a  larger  range  of  incident  angles  where  Eqs 
(64)  and  (78)  do  not  give  the  same  results.  Since  Eq  (78)  is  valid  for 
these  rougher  surfaces  as  well  as  for  the  surfaces  with  small  surface 
slopes,  it  will  be  the  solution  used  throughout  the  remainder  of  this 
thesis . 
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III.  ANALYSIS 


This  chapter  is  devoted  to  numerical  computations  and 
interpretation  of  the  results.  The  sections  in  this  chapter  consist  of 
Part  A: 

1)  A  comparison  of  the  three  pdf’s  presented  in 
Chapter  I  with  respect  to  their  ability  to  represent 
rough  surfaces . 

2)  The  calculation  of  the  average  RCS  per  unit  area  from  a 
surface  represented  by  the  exponential  pdf  (Eq  (3))  along 
with  a  comparison  of  these  results  to  those  for  the  Gaussian 
pdf  and  for  the  exponential-like  pdf. 

Part  B: 

3)  The  computer  generation  of  surfaces  meeting  various 
statistics . 

4)  The  computed  average  RCS  for  the  computer  generated  surfaces. 

Part  A 

1 .  Three  PDF’s 

In  this  section  the  three  pdf’s  presented  in  Chapter  I  will  be 


analyzed  with  respect  to  their  ability  to  represent  rough  surfaces.  .As 
stated  in  Chapter  I,  two  of  the  probability  density  functions,  the 
Gaussian  and  the  exponent ial- type ,  have  been  used  in  the  past  by 


workers  studying  the  electromagnetic  scattering  from  rough  surfaces.  Of 
these  two  pdf’s,  it  will  be  shown  that  only  the  Gaussian  pdf  can 
reasonably  be  expected  to  represent  a  realizable  rough  surface. 

This  author,  along  with  others  [4]  defines  a  realizable 
rough  surface  as  one  for  which  decorrelation  of  the  surface  height 
random  variables  implies  statistical  independence.  Some  may  argue  that 
simply  because  the  variables  do  not  become  independent  does  not  infer 
that  either  variable  influences  the  other  in  a  manner  which  precludes 
the  use  of  the  pdf  from  representing  logical  surfaces.  Rather  than 
contest  this  argument  the  exact  influence  of  one  variable  over  the  other 
will  be  examined  using  conditional  probabilities.  Since  the  bivariate 
Gaussian  pdf  is  a  well  documented  density  which  decorrelates  into  two 
statistically  independent  monovariate  Gaussian  pdf’s,  it  will  be  looked 
at  first  and  used  as  the  basis  of  comparison. 

1 . 1  The  Joint  Gaussian  pdf 

There  are  several  ways  of  determining  the  statistical  independence 
of  random  variables  in  a  multivariate  pdf.  A  sufficient  condition  given 
by  Davenport  is  that  when  the  correlation  coefficient  is  equal  to 
zero,  the  joint  distribution  cam  be  written  as  a  product  of  the  marginal 
distributions  [7:167].  Integrating  the  joint  Gaussian  pdf  with  respect 
to  it  is  found  that  the  marginal  distribution  on  t  is 
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Since  the  joint  distribution  is  symmetric,  t  can  be  substituted  for  t 

i  <1 

in  the  above  equation  to  get  the  marginal  distribution  on  c  .  Setting  p 
equal  to  zero  in  Eq  ( 1 ) 


pG(S'S;0)  =  pG(t2) 


thus  establishing  the  fact  that  decorrelation  implies  statistical 
independence  given  the  joint  Gaussian  pdf. 

Figure  4  illustrates  the  effect  of  the  correlation  coefficient 
upon  the  joint  Gaussian  pdf.  With  P  set  to  zero  ,  Fig.  4a,  one 
simply  gets  the  symmetric  Gaussian  curve  which  states  that  the 


probability  of  <t=  zi  and  <2=  z2  is  dependent  only  upon  the  magnitude  of 
2  2 

{zi  +  z2  ).  When  the  points  are  closer  together  there  is  a  dependence 
between  C  and  C  such  that  the  chance  of  z  being  equal  to  z  is 
increased  while  the  probability  that  zj  will  differ  significantly  from  z2 
is  decreased.  This  trend  is  shown  in  Figs.  4b  and  4c  where  P  is 
equal  to  .5  and  .95  respectively.  With  P= . 95  it  is  seen  that  the 
probability  that  z2is  very  high  as  it  should  be  since  the  points  are 
very  close  to  each  other.  When  P=1  the  joint  density  function  becomes 


6<  Ct-<2)  P  <V 


indicating  that  the  two  random  variables  represent  the  same  point  [4]. 

.Another  way  to  state  statistical  independence  in  the  context  of 
surface  height  random  variables  is  to  say  that  the  probability  of  a 
point  (x,y)  being  at  height  zQ  is  totally  independent  of  the  height  of  a 


point  at  (x  ,y')  as  long  as  the  two  points  are  separated  by  a  great 
distance.  Using  the  definition  of  conditional  probability  [7:160] 
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where  gives  the  probability  of  based  on  a  given  <2  and  p,  a 

necessary  condition  for  statistical  independence  can  then  be  written  as 


p(^  K2;p)  =  P«t) 


(84 


Substituting  the  joint  and  marginal  Gaussian  pdf’s  into  Eq  (83)  gives 
the  conditional  Gaussian  density  function 
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which  is  graphed  out  in  Fig.  5  for  the  same  values  of  P  as  was  the 
joint  pdf. 

With  P-0  the  relationship  of  Eq  (85)  is  illustrated  that  the 
probability  of  '■  is  not  a  function  of  <2>  With  P- .  5  the  center  of  the 
curve  tends  towards  the  <2  line  (it  is  actually  at  the  <  =  .  5t'2 


line)  while  the  variance  of  the  curve  decreases.  With  P-. 95  knowledge 


& 


of  K  almost  completely  determines  ^  as  the  trend  towards 


pVjc^p)  =  6( 


(86) 


becomes  apparent. 

1 .2  Exponential-Type  pdf 

Having  established  some  of  the  desired  features  in  a  joint  pdf  it 
is  time  to  examine  the  exponential -type  pdf  Eq  (2)  which  has  been 
used  by  many  as  an  alternate  to  the  joint  Gaussian  pdf  for  rough  surface 
representation.  Figure  6  graphs  the  exponential -type  pdf  for 
P-0,  .5  and  .95.  There  are  many  similarities  between  these  graphs  and 
the  corresponding  graphs  for  the  joint  Gaussian  pdf .  It  is  these 
similarities  that  may  have  lead  others  to  use  this  pdf  as  a  rough  surface 

MV 

_•  model.  What  is  not  readily  apparent  from  these  graphs  is  that 

the  joint  exponential-type  density  function  is  not  equal  to  the  product 
of  its  marginals  at  P=0  and  and  thus  the  decorrelated  random  variables 
are  not  statistically  independent. 

The  fact  that  the  exponential-type  jpdf  is  not  equal  to  the  product 
of  its  marginals  at  P=0  is  easily  verified  by  checking  at  the  point 
C2=0.  The  marginal  pdf  as  derived  by  Papa  and  Lennon  [15:59]  is 
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where  is  a  first  order  modified  Bessel  function.  Setting  h=l  and 
using  the  value  of  xKi(x)=l  for  x=0  from  [1:379]  one  gets 
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p  (^=0)  = 
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which  does  not  equal  3/(271),  the  value  of  the  jpdf  evaluated  at  zero. 

The  more  critical  question  of  the  relationship  between  C  and  t  is 
answered  in  Fig.  7  which  graphs  the  conditional  pdf’s  for  the 
exponent ial-type  jpdf  at  the  different  values  of  p.  The  fact  that  all 
of  the  graphs  show  a  spreading  of  the  conditional  density  function  on  c 
as  |t2|  gets  large  is  a  curiosity  but  does  not  prohibit  the  jpdf  from 
representing  rough  surfaces.  The  fact  that  the  conditional  density 
function  at  P=0  shows  that  the  density  on  <  is  a  strong  function  of  < 
does  prohibit  its  use. 

The  exponential-type  conditional  pdf  at  P-0  is  odd  indeed  for  a  pdf 
that  is  to  represent  rough  surfaces.  The  increase  in  the  variance  of  C 
as  1 1'2  |  increases  indicates  that  if  the  height  C(x,y)  is  very 
high,  the  probability  that  the  height  at  a  point  (x’,y*)  which  is  very 
far  from  point  (x,y)  is  modified  so  that  very  large  positive  and 
negative  heights  become  more  probable  at  the  point  ( x ’ , y ’ ) .  Meanwhile, 
heights  around  the  mean  become  less  probable.  The  same  behavior  is 
observed  if  the  height  at  point  (x,y)  is  very  large  and  negative. 
Although  this  author  can  offer  no  mathematical  proof  as  to  why  this  is 
impossible  for  physical  surfaces,  it  seems  that  this  is  the  only  logical 
conclusion. 
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Fig.  Ta.  Conditional  Expo-type  pdf,  P=0 
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Fig.  7b.  Conditional  Expo-type  pdf,  P - 
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Fig.  7c.  Conditional  Expo-type  pdf,  P=.95 


1.3  True  Exponential 


The  final  jpdf  to  be  evaluated  is  the  exponential  pdf  introduced  by 
Pyati  ,  Eq  (3),  and  it  is  pictured  in  Fig.  8.  A  noticeable  distinction 
of  this  pdf  is  that  it  does  not  have  the  circular  symmetry  at  p- 0.  This 
is  simply  because  at  P-0  the  exponential  JPDF  reduces  to  the  product  of 
its  marginals 
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which  has  constant  probability  along  the  j t  | + | K2  | =  C  lines,  where  C  is 

.2  .2 

any  constant  instead  of  the  s  +  s  =  C  lines  as  did  the  previous  two  pdf's. 
The  same  trend  as  seen  in  the  other  pdf’s  as  p  approaches  one  is 
apparent  in  Figs.  8b  and  8c.  Figure  8c  is  perhaps  not  as  clean  as 
it  should  be  based  on  the  fact  that  as  P  approaches  1 ,  a  large  number  of 
the  higher  ordered  modified  bessel  functions  are  needed  to  calculate 
accurate  data  points  for  the  exponential  pdf.  Since  each  of  these 
modified  bessel  functions  is  represented  by  a  truncated  infinite  series 
in  the  computer  program  that  calculated  the  surface  graphs,  increased 
error  is  present  as  more  terms  are  needed.  The  trend,  however,  is 
obvious  and  the  fact  that  this  series  tends  to 
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as  P  -»  1  is  verified  when  the  characteristic  function  of  Eq  (3)  is  shown 
to  be  equal  to  1  at  P  =  1  (see  Fig.  12). 
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The  graphs  of  the  conditional  probability  of  £  given  <2  and  P 
(Fig.  9)  show  the  necessary  condition  that  when  P-0  the  probability  of 
<  is  independent  of  <2>  The  above  stated  relationships  as  P  increases 
are  evident  in  Figs.  9b  and  9c  within  the  limits  of  the  numerical 
inaccuracy  stated  in  the  previous  paragraph. 

In  summary,  it  has  been  shown  that  the  true  exponential  JPDF 
satisfies  the  condition  of  decorrelating  into  statistically  independent 
pdf’s  while  the  exponential-type  pdf  has  been  shown  to  be  inappropriate 
for  representing  rough  surfaces. 
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12.  RCS ’ s  Using  Specified  Density  Functions 

In  this  section  the  steps  required  to  calculate  the  average  diffuse 
RCS’s  per  unit  area  are  presented  for  all  three  jpdf’s.  Henceforth  the 
fact  that  only  the  diffusely  scattered  power  is  being  considered  will  be 
assumed  unless  otherwise  stated. 

2.1  The  Characteristic  Functions 

To  calculate  the  average  RCS  per  unit  area  of  a  rough  surface  given 
a  specific  pdf,  one  must  first  calculate  the  characteristic  functions  of 
the  joint  pdf  and  the  one  dimensional  marginal  pdf’s.  The  joint 
characteristic  function  of  a  bivariate  pdf  is  given  by  [7:428] 

*2<vi»v2)  =  <exp(jvt<t+  jv2C2)>  (91) 

Making  the  backscatter  substitution  that  v2  =  -v^:  -v^  this  becomes 

X  (v  ,-v  )  =  <exp( jv  C  -  jv  C  )>  (92) 

2  Z  Z  Z  1  Z  Z 

which,  for  convenience,  will  henceforth  be  referred  to  simply  as  the 
joint  characteristic  function  without  reference  to  the  particular 
arguments.  The  one  dimensional  characteristic  functions  are  given  by 

t  (v  )  =  <exp( jv  <  )>  and  X  (v  )  =  l  ( - v  )  =  <exp(-jv  t  )>  (93) 

1  z  z  lzlz  ^  z 

Before  investigating  the  characteristic  functions  of  specific  pdf’s 
several  properties  of  the  characteristic  function  should  be  noted.  The 
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first  involves  the  case  where  the  two  random  variables  of  a  joint  pdf 


& 


become  completely  correlated.  In  this  case  K  approaches  the  value  of 
<2  and  thus  Eq.(92)  simply  approaches  the  expected  value  of  unity.  The 
second  case  considered  is  when  the  two  random  variable  become 
completely  decorrelated.  If  the  random  variables  become  statistically 
independent  while  decorrelated  the  following  relationship  holds  true 


Lim  t  (v  ,-v  )  =  <exp( jv  K  )><exp(-jv  <  )>  =  r  (v  )x  (v  ) 
P~-*J  i  z  z  z  z  lzlz 


(94) 


If  the  density  functions  are  symmetric  about  zero,  all  of  the  pdf’s  in 
this  thesis  are,  then 
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The  characteristic  functions  of  the  exponential-type  pdf  were  derived  by 
Papa  and  Lennon  and  are  given  as  [15:59] 
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Graphs  of  the  one  and  two  dimensional  characteristic  functions  for  the 
Gaussian  and  exponential -type  pdf’s  are  shown  in  Figs.  10  and  11 
respectively.  The  fact  that  the  exponential-type  does  not  split 
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into  statistically  independent  random  variables  is  obvious  from  the  fact 
that  the  two  dimensional  characteristic  function  does  not  approach  the 
value  of  the  one  dimensional  squared  in  Fig.  11. 


2.1.2  True  Exponential  pdf 

The  one  dimensional  characteristic  function  for  the  true 
exponential  is  easily  calculated  to  be 
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(99) 


Unfortunately,  the  two  dimensional  characteristic  function  is  neither  as 
easy  to  find  nor  is  its  final  form  quite  as  compact.  Because  its 
derivation  is  lengthy  it  has  been  provided  in  a  separate  appendix, 
Appendix  A.  The  solution,  arrived  at  in  collaboration  with  Reif  [19] 
and  with  the  assistance  of  thesis  advisors  Pyati  and  Baker,  is  given  by 
Eq  ( A—  37)  which  contains  two  infinite  series. 

Obviously,  one  would  prefer  to  have  the  solution  given  in  closed 
form.  Appendix  B,  also  derived  with  those  named  above,  gives  a  closed 
form  solution  for  the  first  series  of  the  total  solution.  The  geometric 
series  identity  used  to  get  a  closed  form  solution  for  the  first  series 
cannot  be  applied  to  the  second  series.  The  most  that  could  be  done 
with  the  second  infinite  series  was  to  approximate  it  by  a  finite  series 
and  a  remainder  term.  This  approximation  can  be  found  in  Appendix  C. 

Figure  12  was  generated  by  incorporating  the  contents  of  Appendixes 
A  through  C  into  the  program  MAIN2.BAS  (see  Appendix  D)  which 
calculates  values  of  the  two  dimensional  characteristic  function.  All 
three  characteristic  functions  are  shown  together  in  Fig.  13  using  the 
value  hv  =  10. 
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Correlation  Lengths 
2 

Fig.  12a.  Exponential  Joint  CF  and  *  as  a  Function 
of  Correlation  Lengths  with  h*v2=  2 
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Fig.  12b.  Exponential  Joint  CF  and  ?  as  a  Function 
of  Correlation  Lengths  with  h*v  =  5 


2.2  Performing  the  Final  Integration 


Having  computed  the  characteristic  functions  for  the  pdf’s  all  that 
remains  is  to  perform  the  integration  of  Eq  (78).  Closed  form  solutions 
to  Eq  (78)  have  been  found  in  the  past  for  the  Gaussian  and 
exponential-type  pdf's  where  the  surface  is  assumed  to  be  very  rough 
using  the  method  of  stationary  phase.  Very  rough  in  this  case  means 
that  the  square  of  the  one  dimensional  characteristic  function  is 
essentially  zero  and  the  power  is  diffusely  scattered.  From  reviewing 
Figs.  10  through  13,  this  means  that  the  product  h*v^  should  be  at  least 
5  for  the  Gauss ian  pdf  and  somewhat  larger  than  5  for  the 
exponential- type  pdf.  By  assuming  this  condition,  the  square  of  the  one 
dimensional  characteristic  function  need  not  be  subtracted  from  the  two 
dimensional  characteristic  function  thus  making  the  stationary  phase 
integration  much  easier.  Under  this  condition,  the  diffuse  scattering 
contribution  to  the  average  RCS  per  unit  area  of  a  surface  represented 
by  the  joint  Gaussian  pdf  is  given  by  [20:724] 


G  sec  (©) 

<RCS>d  =  - - -  exp 
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—  tan2(0) 
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( 100) 


2  2  ^ 
where  s  is  the  mean  squared  surface  slope  given  by  4h  /*■  .  ?-  being  the 

correlation  length  of  the  surface.  (A  more  detailed  explanation  of  the 

mean  squared  surface  slope  term  can  be  found  in  Reif  [19]). 
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For  the  exponential-type  pdf,  the  diffuse  power  contribution  to  the 
average  RCS  is  given  by  [20:725] 


<RCS> 
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It  should  be  noted  that  Eqs  (100)  and  (101)  are  valid  only-  for  very- 
large  values  of  h*v^  since  the  method  of  stationary  phase  (a  high 
frequency  asymptotic  technique)  was  used  to  perform  the  integrations. 
This  masks  the  problems  associated  with  the  exponential-type  jpdf.  The 
associated  graphs  of  Eqs  (100)  and  (101)  are  given  in  Figs.  14  and  15 
respectively. 

A  better  physical  understanding  of  Figs.  14  and  15  can  be 
*  obtained  from  considering  a  few  ray  optics  arguments.  As  a  plane  wave 

impinges  normal  to  the  mean  of  a  randomly  rough  surface,  the  rays  will 
be  scattered  in  random  directions  (diffusely).  The  rays  scattered  in 
any  given  direction  can  be  attributed  to  scattering  from  the  specular 
points  for  that  direction.  Increasing  the  roughness  of  the  surface  by- 
increasing  h  does  not  change  the  number  or  position  of  the  specular 
points.  It  does,  however,  decrease  the  principle  radii  of  curvature 
(R  and  R; )  of  those  points.  Since  the  RCS  of  an  object  using  specular 
point  theory  is  given  by  RCS  =  n  RR  [10:119]  this  obviously  will 
decrease  the  average  RCS  of  the  surface.  The  rougher  the  surface  is, 
the  smaller  the  radii  of  curvature  these  specular  points  will  have. 
Therefore,  an  increase  in  roughness  reduces  the  normal  backscatter. 

I 
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.4s  the  angle  of  incidence  is  increased,  the  number  of  specular 
points  is  decreased.  If  the  surface  is  only  slightly  rough,  it  is 
highly  unlikely  that  there  will  be  specular  points  for  angles  of 
incidence  other  than  those  near  zero.  The  rougher  the  surface,  the  more 
likely  it  is  that  there  will  be  specular  points  at  higher  angles  of 
incidence.  {It  should  also  be  remembered  that  the  rougher  the 
surface  is,  the  greater  will  be  the  effects  of  shadowing  which  is  not 
accounted  for  in  these  figures.)  For  more  on  the  specular  theory  of 
rough  surface  scattering  one  can  refer  to  [4],  Reif  [19]  is  a  good 
source  for  the  shadowing  theory. 

Before  taking  on  the  task  of  performing  the  numerical  integration 
of  Eq  (78)  with  data  points  of  the  exponential  characteristic  functions, 
the  behavior  of  the  elements  of  the  integrand  should  be  analyzed  with 
respect  to  their  possible  effects  on  the  integration.  Keeping  in  mind 
what  must  be  done  in  a  numerical  integration,  it  is  seen  that  to 
reproduce  solutions  such  as  those  in  Figs.  14  and  15,  one  must 
successively  solve  the  integral  for  discrete  values  of  the  mean  surface 
slope  and  at  discrete  angles  of  incidence.  The  variables  that  change  as 
these  parameters  change  are  h,  v  and  while  the  correlation  length, 
t,  is  kept  constant.  (The  correlation  length  is  set  to  1  throughout 
this  thesis.)  Two  of  these  variables,  h  and  v  are  arguments  of  the 
characteristic  functions.  Because  of  the  length  of  time  it  takes  to 
calculate  tiata  points  for  the  exponential  characteristic  function  (  To 
generate  the  5001  data  points  that  were  eventually  used  in  the 
integration  it  took  56  hours.)  one  must  find  a  way  of  using  only  one  set 
of  data  points  from  the  characteristic  function  and  scaling  the  other 


parameters  accordingly.  Fortunately,  h  and  v  only  appear  as  the 
dimensionless  parameter  h*vz  within  the  characteristic  function. 
Therefore,  the  goal  can  be  accomplished  by  letting  h*v  equal  some 
constant . 

To  keep  h^v^  constant  for  any  given  mean  surface  slope  as  the  angle 
of  incidence  changes,  the  following  must  hold: 


h*v  =  2k  cos(0)  =  C 

z 


( 102) 


k  = 


C 

2h  cos(@) 


(  103) 


where  C  is  the  yet  to  be  determined  constant.  Fortunately,  as  can  be 
seen  by  reviewing  Eqs  (100)  and  (101),  the  <RCS>D  is  not  a  function  of  k 
as  long  as  k  is  large  enough  so  that  the  one  dimensional  characteristic 
functions  can  be  ignored.  From  Fig.  13  it  would  appear  that  letting 
C  =  10  is  sufficiently  large,  therefore  10  will  be  the  value  used.  The 
use  of  C  :  10  as  a  high  enough  number  will  also  be  confirmed  by  the 
accuracy  with  which  the  numerical  integration  of  the  <RCS>  integral 
reproduces  the  stationary  phase  result  for  the  Gaussian  surface. 

It  is  seen  from  Eq  (103)  that  as  the  angle  of  incidence  increases, 
the  value  of  k  to  be  used  in  the  integration  increases.  This  causes  the 
constant  portion  of  the  Bessel  function  argument  in  Eq  (78)  to  be 
increased.  This  in  turn  increases  the  frequency  of  the  Bessel  function. 
Figure  16  illustrates  this  relationship.  The  fact  that  it  is  the 
increased  oancellation  between  the  positive  and  negative  portions  of  the 
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high  frequency  Bessel  function  that  cause  the  RCS  to  go  to  zero 
emphasizes  the  necessity  to  use  enough  data  points  in  the  integration. 
When  the  angle  of  incidence  is  held  constant  and  the  mean  surface  slope 
is  increased,  the  result  is  a  decrease  in  the  wave  number.  The  effects 
of  this  relationship  are  shown  in  Fig.  17. 

Since  the  characteristic  functions  for  the  Gaussian  and  exponential 
pdf’s  are  not  dramatically  different,  one  might  assume  that  the  same 
number  of  data  points  required  to  integrate  Eq  (78)  using  the  Gaussian 
CF  would  be  adequate  to  integrate  the  exponential.  By  incrementally 
increasing  the  number  of  data  points  used  to  do  such  an  integration,  it 
was  found  that  with  h*vz=  10  that  5001  data  points  evenly  spaced 
from  0  to  1  correlation  length  (  that  is,  .0002  correlation  lengths 
apart  )  produced  adequate  results.  These  results  are  shown  in  Fig.  18. 
Overlaying  this  graph  and  the  graph  of  the  stationary  phase  solution  of 
the  integral,  it  is  found  that  the  values  are  exactly  correct  until  they 
get  down  into  the  -20  dBsm  region  where  they  obviously  deviate  from  the 
stationary  phase  solution. 
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Creating  the  5001  data  points  for  the  integration  proved 

to  be  no  trivial  task.  As  stated  earlier,  every  attempt  was  made  to 

minimize  the  number  of  calculations  required  to  produce  this  exponential 

CF  data  (see  Appendixes  A-C).  Finally,  the  program  MAIN2.BAS  was 

written  to  crunch  the  numbers.  It  was  soon  realized  that  running  this 

program  on  the  basic  IBM  PC,  even  in  compiled  TurboBASIC,  would  place 

the  run  time  of  the  program  dangerously  close  to  the  mean  time  between 

failures  of  the  computer.  The  computer  was  therefore  enhanced  with  a 

faster  CPU  and  a  numerical  coprocessor.  This  increased  the  speed  of  the 

calculations  by  a  factor  of  8  allowing  the  data  file  to  be  built  in  only 

56  hours  of  continuous  run  time.  After  completing  the  data  file, 

program  RCS2.BAS  was  run  to  perform  the  integration  using  a  Simpson’s 

integration  scheme.  The  results  of  these  integrations  were  then 
2 

multiplied  by  sec  (0)  to  produce  Fig.  19. 

Figure  19  gives  the  average  RCS  per  unit  area  of  a  very  rough 

surface  represented  by  the  true  exponential  pdf.  It  is  interesting  to 

note  that  it  does  indeed  resemble  the  graph  of  the  exponential-type 

results  in  shape.  This  is  the  only  similarity  though,  as  overlaying  the 

graphs  clearly  shows  that  the  values  are  not  the  same.  The  RCS  of  the 

true  exponential  surface  is  lower  at  normal  incidence  and  falls  off  much 

slower  than  that  of  either  the  exponential-like  or  Gaussian  surfaces.. 

This  indicates  that  the  true  exponential  surface  is  a  rougher  surface 

2 

given  the  same  surface  height  variance  (h  ).  In  the  next  section, 
computer  generated  surfaces  will  show  this  to  be  true  but  for  now  this 
argument  will  be  supported  by  simply  looking  at  a  comparison  of  the 


3 .  Summary 


In  this  section  some  of  the  properties  of  the  Gaussian, 
exponential-type  and  the  exponential  joint  probability  density  functions 
were  investigated  along  with  some  of  the  related  properties  of  their 
characteristic  functions.  It  was  shown  that  upon  decorrelation  the 
Gaussian  and  exponential  random  variables  become  statistically 
independent.  It  was  also  shown  that  the  exponential-type  jpdf  lacks 
this  property  and  is  therefore  not  suitable  as  a  rough  surface  model. 

In  support  of  the  <RCS>  results  found  for  the  different  jpdf 
surface  representations  one  can  refer  to  Figures  4,  6  and  8  all  of  which 
were  generated  with  the  variance  set  equal  to  one.  Although  the  scale 
of  the  graphs  are  all  different  (this  was  done  to  maintain  constant  size 
for  the  graphs)  it  is  obvious  that  the  exponential-like  and  the  Gaussian 
pdf’s  go  to  zero  much  faster  than  does  the  exponential.  This  indicates 
that  the  exponential  surface  will  have  higher  peaks  on  its  surface  than 
will  the  other  surfaces.  This  in  turn  will  give  the  exponential 
surfaces  smaller  specular  points  distributed  over  a  larger  range  of 
angles.  Referring  back  to  the  specular  point  arguments  for  the  way  the 
RCS  of  a  rough  surface  acts,  this  supports  the  findings  presented  in 
Figs.  14,  15  and  19. 
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Part  B 

Part  A  of  this  section  dealt  with  mathematically  predicting  the 
average  diffuse  RCS  of  a  surface  given  that  the  surface  was  represented 
by  specific  probability  density  functions.  The  only  way  to  check  such 
calculations  would  seem  to  be  to  check  them  against  RCS  measurements 
made  on  real  surfaces  that  obey  such  statistics.  Most  surfaces  in 
nature,  however,  have  many  different  degrees  and  types  of  roughness 
superimposed  upon  one  another.  Therefore,  surfaces  that  conform  to  a 
single  statistical  representation  are  not  easily  found.  (Two  references 
that  address  this  topic  are  [13]  and  [18].)  An  alternative  given  in 
this  section  to  real  RCS  measurements  is  computer  predicted  RCS’s  of 


computer  generated  rough  surfaces  that  meet  the  proper  statistics. 

In  this  section,  random  rough  surfaces  that  obey  Gaussian  and 
exponential  statistics  will  be  generated  by  computer.  The  RCS  of  these 
surfaces  will  then  be  calculated  for  different  combinations  of  mean 
surface  slopes  and  h*v^.  These  results  will  then  be  graphed  and 
compared  with  the  predicted  results  from  Part  A. 

1 .  The  Generation  of  Random  Rough  Surfaces  ( Background) 

The  technique  for  producing  random  rough  surfaces  is  borrowed  from 
Karnly  and  Jensen  [9]  where  it  was  used  to  create  spatially  correlated 
infrared  noise.  Harnly  and  Jensen,  developed  this  technique  under  the 
guidance  of  AFIT  faculty  member  Dr.  Peter  Maybeck.  The  only 
modification  to  the  technique  is  one  of  interpretation.  That  is, 
instead  of  creating  spatially  correlated  noise,  one  is  now  finding 


spatially  correlated  surface  heights. 
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The  correlation  matrix  is  simply  an  N  x  N  matrix,  each  element  of 
which  represents  the  correlation  coefficient  between  a  point  (j)  and  a 
point  (k)  where  (j)  and  (k)  index  from  1  to  N2  [7:366].  The  correlation 
matrix  is  given  by 


For  the  purposes  of  this  thesis,  a  Gaussian  correlation  coefficient  is 
used.  Therefore 

P  .=  exp  (-d2  /f 2  )  (106) 

j.k  j,x 

where  d^  k  is  the  distance  between  points  j  and  k.  Obviously,  it  then 

follows  that  P .  .  =  1  and  P  ,  =  P,  ,  .  Also,  as  the  distance  between 

j.j  j.x  X.j 

two  points  gets  large,  the  correlation  coefficient  will  go  to  zero  and 
the  points  will  be  statistically  independent.  The  extreme  example  of 
this  would  be  a  grid  for  which  even  the  adjacent  points  were  very  far 
apart.  Under  this  condition  the  correlation  matrix  would  simply  be  the 
identity  matrix  which  has  off  diagonal  terms  of  zero  and  diagonal  terms 
equal  to  one.  In  such  a  case  Eq  (104)  simply  gives  Z(t  )=W(t  ). 


The  elements  of  the  Cholesky  square  root  matrix  are  generated 
sequentially,  row  by  row,  from  the  recursion:  for  i = 1 , 2 , . . .  X2  compute 


[12:371] 
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where  the  a’s  are  the  elements  of  the  Cholesky  square  root  matrix  and  the 

P’s  are  the  correlation  matrix  elements. 

—  2 

As  stated,  the  vector  W(t  )  is  simply  N  uncorrelated  surface 
heights  produced  by  any  valid  one  dimensional  pdf.  This  is  extremely 
convenient  since  one  dimensional  pdf’s  are  easy  to  come  by.  Surfaces 
can  then  be  generated  which  follow  Gaussian  statistics,  exponential 
statistics  or  even  uniform  or  Raleigh  statistics  if  one  so  desired. 

The  program  that  will  be  used  to  calculate  the  RCS  of  the  generated 
surfaces  is  written  in  Fortran  and  resides  on  a  VAX-11/785  which  is 
owned  by  AFWAL.  Therefore,  the  programs  that  generate  the  surfaces  were 
written  in  Fortran  on  the  same  computer.  The  following  programs  were 
written  for  this  purpose 

CORCOEFF :  Given  the  horizontal  distance  between 
surface  points  and  the  number  of  points  in  a  row, 


this  program  calculates  the  total  distance  between 
each  point  and  creates  a  file  containing  the 
correlation  matrix. 


CHOLESKY:  Given  a  correlation  matrix  file,  this 
program  creates  the  lower  triangular  Cholesky  square 
root  of  it  aind  plaices  it  in  a  file. 

GESURF:  Given  a  Cholesky  square  root  file,  this 
program  creates  the  W  ( b  )  vector  and  then  the  Z  ( b  ) 
vector  which  represents  the  rough  surface.  The 
surfaces  generated  by  this  program  have  a  surface 
height  variance  of  one.  The  variance  of  these 
surfaces  cam  be  changed  by  simply  multiplying  the 
existing  heights  by  the  square  root  of  the  variance,  h. 

Each  of  these  programs  is  listed  in  the  Appendix  D. 

The  diagonal  terms  in  the  Cholesky  decomposition  algorithm  get 
increasingly  smaller  as  the  size  of  the  correlation  matrix  increases  and 
the  horizontal  distance  between  points  increases.  Roundoff  error  then 
comes  into  play  and  can  cause  the  term  within  the  square  root  sign  on 
Eq  (108)  to  become  negative  which  should  not  happen  for  the  positive 
definite  correlation  matrices  presented  in  this  thesis.  This  limits  the 
the  size  of  the  surfaces  that  can  be  generated.  Trial  and  error  showed 
that  square  surfaces  could  be  generated  that  were  a  total  of  1.2 
correlation  lengths  long  and  had  15  rows  and  columns  of  surface  heights 
each  separated  by  .3  correlation  lengths.  The  generated  surfaces  are 
shown  in  Figs.  El  through  E20  of  Appendix  E  along  with  the  RCS 
information  for  each  surface. 
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2 .  The  RCS  Program 

The  program  used  to  calculate  the  RCS  of  the  computer  generated 
surfaces,  the  Radar  Cross  Section-Basic  Scattering  Code  (RCS-BSC),  was 
provided  by  AFWAL/AAWP-3 .  The  code  itself  was  written  at  Ohio  State 
University's  Electroscience  Laboratory  under  the  supervision  of  Ronald 
Marhefka  [11].  Although  the  RCS-BSC  has  capabilities  well  beyond  that 
of  just  providing  physical  optics  predictions  for  perfectly  conducting 
bodies,  only  the  option  which  makes  these  assumptions  was  used  since  the 
calculations  in  Part  A  are  based  upon  these  assumptions. 

Two  minor  problems  had  to  be  overcome  before  the  RCS-BSC  could  be 
run  against  the  generated  surfaces.  First,  the  RCS-BSC  could  only 
handle  objects  represented  by  a  maximum  of  50  plates  with  a  maximum  of 
40  sides  per  plate.  In  addition,  the  code  does  not  allow  plates  to  be 
concave.  That  is,  all  of  the  corners  of  the  plates  must  be  in  the  same 
plane.  The  surfaces  generated  are  made  up  of  196  plates,  each  with  4 
sides,  none  of  which  are  planar.  Therefore  RCS-BSC  had  to  be  modified 
to  allow  objects  to  be  represented  by  more  plates.  This  was  simply  a 
matter  of  plodding  through  the  program  and  changing  the  appropriate 
dimension  statements.  Since  no  more  them  foirr  sides  per  plate  were 
needed,  it  was  possible  to  increase  the  number  of  plates  the  program 
could  handle  while  decreasing  its  memory  requirement  by  redimensioning 
the  plate  geometry  array  to  400  -  4.  The  problem  of  concaveness  was 
solved  by  representing  each  square  plate  as  two  triangles  since 
triangles  are  planar  by  definition.  This  would  have  allowed  the 
red Lmensioning  of  the  arrays  to  allow  only  three  sides  instead  of  four 
but  this  was  not  done  since  computer  memory  availability  was  not  a 


problem  and  since  the  redimensioning  to  allow  up  to  400  plates  had 
already  been  accomplished.  Since  each  of  the  196  squares  was  now- 
represented  by  2  triangles,  382  of  the  400  allowable  plates  were  used. 

To  finally  calculate  the  RCS  of  the  surfaces,  files  needed  to  be 
generated  that  set  the  proper  control  parameters  for  RCS-BSC  and  placed 
the  surface  height  data  in  the  proper  format.  Program  BUILD. FOR 
accomplished  this  chore.  In  addition,  since  a  single  run  of  RCS-BSC 
against  a  surface  took  27  CPU  minutes  on  the  VAX,  batch  files  had  to  be 
written  to  run  the  code  after  hours. 

As  mentioned,  the  run  time  of  the  program  was  a  concern.  Even 
limiting  the  RCS  calculations  to  5  degree  increments  between  0  and  90° 
the  calculation  of  a  single  pattern  took  27  CPU  minutes.  Therefore, 
only  a  limited  number  of  surface  RCS  calculations  could  be  run.  The 
final  series  of  runs  consisted  of  10  Gaussian  surfaces  and  10 
exponential  surfaces.  The  heights  of  these  surfaces  was  scaled  so  that 
each  had  mean  surface  slopes  of  5°and  25°.  The  RCS  of  each  of  these 
surfaces  was  in  turn  calculated  at  5.46  Ghz  and  17  Ghz  so  that  the 
pat-tern  could  be  compared  as  a  function  of  the  wavenumber  k.  The  RCS’s 
calculated  for  each  surface  are  given  in  the  tables  associated  with  each 
surface  in  Appendix  E. 

Having  completed  all  the  runs,  the  RCS’s  were  averaged  at  the 
discrete  angles  of  incidence.  To  average  the  RCS’s  of  the  surfaces, 
each  RCS  was  first  converted  to  square  meters.  The  average  of  the  RCS's 
was  then  calculated  in  square  meters  and  then  transformed  into  dBsm. 
Tables  IT  and  III  give  the  results  of  this  process  along  with  the 
physical  optics  solution  for  the  RCS  of  a  flat  plate  at  each  point. 


These  results  are  also  graphed  in  Figs.  21  through  24.  The  RCS’s  of  the 
flat  plates  are  provided  so  that  one  may  see  that  the  diffuse  scattering 
does  indeed  fill  the  nulls  of  the  coherent  pattern  and  decreases  the  the 
peaks  of  the  pattern. 

The  results  of  this  experiment  show  excellent  agreement  with  the 
predicted  results.  The  selected  combinations  of  mean  surface  slopes  and 
frequency  assured  that  hk  >  5  for  all  measurements.  Therefore,  as 
stated  in  Part  A,  the  scattered  power  should  be  primarily  diffuse  and 
independent  of  frequency.  Comparison  of  the  results  shows  that 
frequency  independence.  In  addition,  comparison  of  Figs.  (21)  through 
(24)  with  Figs.  (18)  and  (19)  show  excellent  agreement  between  the 
behavior  of  the  average  RCS  as  a  function  of  mean  surface  slope  and 
incident  angle  as  well  as  good  agreement  between  the  actual  RCS  values 
and  the  predicted  RCS  values  at  small  angles  of  incidence. 
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Fig.  21.  <RCS>/Area  for  Computer  Generated  Gaussian 
Surfaces  at  5.46  Ghz 
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Fig.  22.  <RCS>/Area  for  Computer  Generated  Gaussian 


Surfaces  at  17.0  Ghz 
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Fig.  23.  <RCS>/Area  for  Computer  Generated  Exponential 
Surfaces  at  5.46  Ghz 
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IV.  Conclusion  and  Recommendations 


1 .  Summary  of  Work  Completed 

In  this  thesis  it  was  postulated  that  any  joint  probability 
density  function  used  to  model  statistically  rough  surfaces  should 
decorrelate  into  statistically  independent  random  surface  height 
variables.  A  joint  exponential  pdf  which  met  this  requirement  was  then 
introduced  and  the  diffuse  scattering  properties  of  surfaces  modeled  by 
it  were  derived  using  physical  optics  and  infinite  conductivity 
approximations.  This  represents  the  first  time  the  diffuse  scattering 
properties  of  a  true  joint  exponential  pdf  have  been  calculated.  The 
resulting  scattering  properties  of  the  joint  exponential  pdf  surfaces 
were  then  compared  to  the  scattering  properties  of  the  joint  Gaussian 
pdf  and  a  joint  exponential  pdf  and  were  found  to  be  significantly 
different. 

In  addition,  a  technique  was  introduced  whereby  one  can  generate 
rough  surfaces  that  satisfy  virtually  any  kind  of  statistics.  Then,  10 
surfaces  each  that  satisfied  Gaussian  and  exponential  statistics  with  mean 
surface  slopes  of  5°  and  25°  were  generated  and  the  high  frequency  RCS’s 
of  these  surfaces  was  calculated.  The  resulting  average  RCS’s  per  unit 
area  for  the  different  types  of  surfaces  were  shown  to  be  in  in  good 
agreement  with  the  RCS’s  predicted  by  the  more  traditional  means  but  to 
deviate  from  the  traditional  predictions  off  normal  angles  of  incidence. 


2.  Recommenda t i ons  for  Further  Stucb 
I  .  '  '  . 

Several  follow  on  studies  that  could  be  conducted  are  enumerated 

below. 

(1)  The  rough  surfaces  studied  in  this  thesis  could  be  expanded  to 
include  rough  surfaces  that  have  finite  conductivity.  The 
calculations  could  be  done  using  Eq  (79)  with  the  scattering 
matrix  elements  defined  for  the  particular  surface  under 
study.  Introduction  of  finite  conductivities  to  the  surfaces 
generated  by  computer  in  Chapter  III,  Part  B  could  be  easily 
accomplished. 

(2)  The  effects  of  different  correlation  coefficients  upon  the 
diffuse  scattering  properties  of  the  joint  exponential  pdf 
could  be  studied.  This  would  require  only  slight  modification 
to  the  integrand  of  Eq  (78)  since  the  characteristic  functions 
are  a  function  of  the  correlation  coefficient.  It  would  also 
be  possible  to  replace  the  Gaussian  correlation  matrix  used  in 
Chapter  III,  Part  B  with  some  other  type  of  correlation 
matrix.  This  was  tried  for  an  exponential  correlation 
function  with  no  apparent  problems. 

(3)  With  the  introduction  of  the  rough  surface  generating 
technique,  surfaces  that  meet  other  type  of  statistics  could 
be  generated  and  the  RCS’s  calculated  for  comparison  to  the 
ones  studied  in  this  thesis. 

(4)  Follow  on  work  could  be  done  to  show  the  validity  of  modeling 
rough  surfaces  with  joint  probability  density  functions.  In 


this  regard,  the  surface  heights  of  actual  terrain  and  other 


rough  surfaces  of  interest  could  be  measured  and  statistically 
analyzed.  Also,  actual  Laser  and  radar  cross  section 
measurements  could  be  made  and  studied  for  comparison  to  the 
predicted  results. 

(5)  Surfaces  appearing  in  nature  usually  have  more  than  one  type 
of  surface  roughness.  For  instance,  terrain  might  have  large 
scale  roughness  that  is  properly  modeled  with  a  certain  joint 
probability  density  function  while  its  smaller  scale  roughness 
might  be  more  properly  modeled  by  another  type  of  jpdf.  On 
top  of  the  terrain  one  would  possibly  find  vegetation  which 
might  best  be  modeled  by  randomly  oriented  cylinders. 
Therefore,  a  study  that  showed  the  construction  of  an  entire 
surface  model  might  be  appropriate. 

(6)  Finally,  a  study  that  looked  at  the  total  bistatic  scattering 
of  different  types  of  rough  surfaces  could  be  done.  This 
study  should  include  an  integration  of  the  total  power 
scattered  diffusely  for  different  scales  of  roughness.  This 
power  could  then  be  subtracted  from  the  total  power  incident 
upon  the  surface.  The  remaining  power  would  be  the  coherently 
scattered  power.  The  total  scattered  power  would  be 
represented  by  the  addition  of  the  diffuse  and  coherent 


scattered  fields. 


The  purpose  of  this  appendix  is  to  show  how  the  characteristic 

E 

function  X AK  ,K2',P)  for  Pyati’s  two  dimensional  exponential  pdf, 

Eq  (3),  was  derived.  The  characteristic  function  ~l2  ( x\  >  v2  )  °f  a 
two  dimensional  pdf  P  I'C  )  is  defined  as  [6:189] 

i  X  Z 

X)  00 

VVV  =  f  1  P2(W  eXPC  jv1<1  +  jv2g  d^d<2  (A- 

-0D  -0D 

To  find  the  characteristic  function  of  Dr  Pyati’s  two  dimensional 
exponential  pdf,  replace  P, (<.,<,)  in  Eq  (A-l)  with  Eq  (3)  and  solve. 
The  integrals  that  need  to  be  solved  to  obtain  the  characteristic 
function  are  shown  in  Eq  (A-2). 


In  Eq  (A-2),  ^  and  <2  symbolize  two  random  surface  heights  which 
are  functions  of  x  and  y.  <  and  <2  are  separated  by  a  distance  r  and 
where  p  •  p(r)  is  the  surface-height  correlation  function.  I  is  a 
modified  Bessel  function  of  the  first  kind  of  order  a.  Also,  sgn(x)  is 
the  signum  function  of  x. 
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At  first  glance,  it  appears  impossible  to  obtain  a  closed  form 
solution  for  Eq  (A-2).  However,  after  much  trial  and  error,  a  valid 
solution  was  obtained. 

The  first  step  in  solving  Eq  (A-2)  is  to  remove  the  dependence  of 
the  integral  upon  the  absolute  value  functions  and  the  signum  functions 
This  can  be  done  by  using  the  following  relationships. 


I <,.,!=  -*%2  for  <x.a  <  0 


'  ,,2  f°r  0  <  <  * 


sgn(<  )  =  -1  for  ->»  <  <  <  0 

1,2 


=  +1  for  0  <  <  < 

1  •  " 


Using  equations  (A-3)  and  (A-4)  the  above  Fourier  transform 
integral,  Eq  (A-2),  can  be  written  as  the  sum  of  four  integrals. 


VW  =  t  +  :2  +  *  h 


The  values  of  I  through  I  are  given  on  the  next  two  pages.  Also, 

1  4 


a  substitution  of  c  =  h(l-P  )  is  used  for  brevitv. 
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These  integrals  need  to  be  arranged  into  a  form  that  can  be 


integrated  or  for  which  an  integral  is  known.  One  form  that  works  was 
(see  Appendix  A.l)  is  given  by 


p-l  -ax  , 

x  e  cos(mx)  dx  = 


T ( p )  cos ( p@ ) 


,  2  2  d/2 

(a  +  m  ) 


(A- 10 ) 


p-l  -ax 

x  e  sm(mx)  dx  = 


T ( p )  s in ( p£ ) 

,  2  2  d/2 

(a  +  m  ) 


( A- 1 1 ) 


0  =  tan 


-i  m 


2  2.1/2 
(a  +  m  ) 


( A- 1 2 ) 


Several  manipulations  must  be  performed  to  get  integrals  one 
through  four  into  the  above  form.  The  first  step  is  to  make  all  the 
limits  of  integration  to  be  from  zero  to  infinity.  This  is  done  in  two 
steps,  first  making  all  integrals  start  at  zero,  then  changing  any 
negative  infinities  to  positive  by  changing  the  sign  of  the  appropriate 
variable.  These  changes  are  shown  on  the  next  four  pages. 


These  integrals  can  be  simplified  by  noting'  that  they  are  really 
the  sum  of  two  integrands.  Thus,  if  the  simplifications  listed  below 


MimflMimM»W«lJWIUUUUUmWHJMYWW'WtJ»tfwiWlHnrwimgnxKJ’mjMi  r*  rrrr.-arsmmnm  wrm.-w»nnwmi  ni 

Equation  (A-26)  with  the  use  of  a  the  trigonometric  identity  for 
the  addition  and  subtraction  of  angles,  which  states  cos(Ux.J)  =  cos« 
cosd  ?  sina  sin  3,  can  be  written  as 


E 

^  = 


2  /  /  D  (A  +  B)  j  cos  cos  (v  <■'  ) 

0  0  1 

-  sin^^J  sin(v2<2]  } 

x>  oo 

+  2/  f  D  (A  -  B)  -j  cos(vi<i3  cos(v2<o3 

sinCv^J  sin(v2<2]  j 


+  sm 


(A-27) 


Expansion  of  Eq  (A-27)  simplifies  to 


od  oo 


o  o 


X  =  4  /  /  D  (A  c  os(v  <  )cos(v  <  )  ]  d£d<, 

2  _  _  11  2  2  12 

00  00 

-  4  /  /  D  (  B  sin(v  <  )sin(v2<  )  ]  dC  d*.'. 


( A-28 ) 


Now  the  following  substitutions  <  =  C  / c  and  >7  =  <  /c  are  made 

1  2 

and  the  Bessel  function  is  replaced  with  its  series  equivalent  Eq  { A— 29 ) 
found  in  [  1 : 375 ] . 


Vz) 


(z2/ok 


k  !  r ( l +k+  1  ) 


k  »  0 


( A- 29 ) 
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When  «-  and  k  are  integers,  r  (t*+k+l)  =  (t+k)!  [1:255]  and  Eq  { A— 29 ) 

reduces  to  Eq  { A— 30 ) . 


Vz)  = 


C  z  2  /  4  )  k 


k!  (v+k)! 


( A-30 ) 


Eq  (A-30)  and  the  <?and  h  substitution  are  used  to  transform  the 
variables  A,  B,  and  D,  given  in  Eq  (A-21)  through  (A-23),  into  the 
following  equations. 


I  C2P 

O 


Cp2<h)k 


( A-31 ) 


1 2n+ 1  ^  3 


(2n  +  I) 


XX  2  n  ♦  1 

■  ■■  \  f*  "  ^  i  /  2  1  /  2  ^ 

v  \  Cp  <  ~n  ) 


(p  th) 


(2n  +  1)  k!  (2n+l+k) 


{ A— 32 ) 


k  »  0  n  ■  0 


exp  C~(  -n) 


( A-33 ) 


K 

Kv 


I 

I 


Substituting  Eq  (A-31)  through  (A-33)  into  Eq  ( A— 28 )  allows 


Eq  (A-34)  now  has  its  integrals  in  the  form  of  equations  ( A— 1 0 )  and 


< A— 11)  for  which  solutions  are  known.  Therefore,  the  characteristic 
£ 

function  t  for  Pyati’s  two-  dimensional  exponential  pdf  has  the 
following  solution. 


=  o2k  r2  (k+1)  cos((k+l)@i]  cos((k+l)^2) 

~  /  .  ,  2  T  "  2*.  c  k  + 1 )  /  2  7"  "  TiTTkTTTTT 

h  l_ ,  k!  [1  +  (cvj  J  [1  +  <cv2)  J 


2n  +  1  +  2  k 


C2n  +  1)  k!  ( 2n+ 1 +k ) ! 


n  »  0  k  »  0 


r2(n+k  +  l.o)  sin  [(  n  +  k+  1 . 5  )  9  ^ )  s  i  n  [(  n  +  k+  1  .  5  )  9  ] 


Ci  +  ( cvx )  3 


2-.<ntk  +  l  .  5)  /  2 


(l  +  (cv  )  ] 


2  \  (  n  +  k  t 1  .  5  >  /  2 


( A-35 ) 


Where  9 ^  and  9 2  are  given  as 


f  *• 

\.2  =  t^1  [  CV1.2  , 


,  2  2  .1/2 
{  1  +  C  V  t  ,  ) 


( 1-P  )hv. 


2  2  1/2 
(  1  +  C  V  ) 

1.2 


;  A-36 ) 


To  obtain  the  variable  0  ,  just  replace  v  with  v  in  Eq  ( B— 2 )  . 

Then,  using  Euler's  identity  for  a  cosine,  the  cosine  term,  cos  £( k+1 )©  ] 
is  written  as  Eq  (B-3).  The  same  is  done  for  the  cos((k+l)©23  term  by 
replacing  0  with  &2  . 


cos( (k+1 )  = 


exp(j  (k+1  )6>i]  +  exp  f -j  (  k+ 1  )  ©1  ] 


Also,  to  reduce  Eq  ( B— 1 )  to  a  workable  form,  the  following  substitutions 
will  be  used. 


A  =  fl  +  (WjhvJ2  11/2 


f  -\ 

I  =  1  +  ((l-p2)hv2)2 


The  use  of  equations  (B-3)  through  ( B— 5 )  allows  equation  ( B— 1 )  to 


be  written  as 


(  1-P 


k  k 

A  B 


exp  (j  (k+1  +  exp  ( -j  (k+1  )©i 


o 

exp(j(k+l  )©  ]  +  exp(-j(k+l)@2]  J- 


t J 


Expanding  Eq  ( B— 6 )  by  multiplying  terms  to  remove  the  brackets 


yields 


U -P  ) 


"v~7~  1  exp[j(<9 +0  )]  exp(jk(0 +0  )) 
A  B  12  11 


+  exp[j(0i-02  )]  exp  Cjk(^-^2)) 

+  exp[-j(01-02  )]  exp(-jk(0i-02  )] 

+  exp  (-  j  ( 0^  +0  2  ) }  exp  (- jk  ( Qy  +02  ) } 


Using  the  geometric  series  identity 


•  T*  -  - 


(B-7) 


which  is  valid 


when  |  X  |  <  1,  [MH:107]  equation  (B-7)  can  be  written  as  Eq  (B-8). 

Equation  (B-8)  is  now  in  a  close  form  since  the  summation  sign  has  been 
removed . 


(  1-P  ) 


exp(j(0i+02  )} 

ACB  -  p2exp  (j  (@i  +@2  )} 


exp  (j  (el~92  >]  exp[-j(0  -0  )] 

'  +  '  '  ' 

•  p2exp  [j  ( &l -&2  )]  AB  -  p2exp  (- j  (  -02  )} 


exp  (-j  (  0t+02  )]  | 

AE  -  p2exp(-j(0i+02)J  J 


(B-8) 


Manipulating  Eq  ( B— 9  >  by  multiplying  each  term’s  numerator  and 


In  Eq  (B-12)  and  ( B— 13),  the  trigonometric  identity  for  the  cosine 
of  the  sum  of  two  angles  is  used  to  expand  J-!  and  F  in  terms  of  sine  and 
cosine  of  #.  Also,  Eq  { B— 2 ) ,  Eq  < B— 4 )  and  Eq  ( B— 5 )  are  used  to  evaluate 


the  sine  and  cosine  of  #. 


iM  =  cos(#  +#  )  =  cos#  cos#  -  sin#  sin# 

1  2  12  12 


F  =  cost#  -#  )  =  cos#  cos#  +  sin#  sin# 
12  12  12 


cos (cos1 ( 1/A ) 


cos (cos  ( 1/E  ) ) 


|  (1-0  !hvi 

sin  sin  - 

A 


j  j  (  1  -P  )  hv2 

sin#  =  sin  sin  - 

1  E 


(1-0  )hv( 


(B-12) 


( b- i 3 : 


( B- 14  ) 


( B- 1 5 ) 


{ B— 1 6  ) 


(1-0  ) hv. 


Equations  (B-12)  through  (B-17)  are  used  to  write  new  equation  for 
and  7  .  These  terms  are  now  void  of  sine  and  cosine  terms. 


1  -  (l-02)h%i  (  1  -  O 2  )  h  x-  ^  1  -  (l-02)2h2viv2 


1  +  (l-02)h\'i  (  1  - 0 2  )  h v ^  1  +  (  1-02  )  2h  2vj 


(B-18) 


(B-19) 


VV. 


I 


# 


3 


.Substituting  equations  (B-18)  and  ( B— 19}  into  Eq  ( B— 1 1 )  yeilds  Eq 


2^2 


(B-20).  The  values  for  A  E  are  given  in  Eq  ( B— 2 1 ) 


V 


( i -p* )  r  i  -  < i-p2  > 2 h 2 v 1 v 2  -  p2 

A2B2  -  2P2  (l-(  1-p2  )  2h2vt  vj  +  P4 


2  2  2  2 
1  +  (1-P  )  h  vjv2  -P 

,  2  2  2  c  .  2  2,2  ^  4 

A  E  -  2P  [l  +  (  1-P  )  h  v  t  v  2  J  +  p 


(B-20) 


with 


A2E2  =  [l+  (  ( \-pL  )hv^  I*1]  (l  +  (  ( \-P  )h\2  j"1] 


( B— 2 1 ) 


Combining  Eq  (B-20)  and  (B-21)  and  factoring  like  terms  reduces  the 
final  equation  to  a  close  form  solution  for  the  first  series  of  Eq 
(B-l ) . 


1 

v- 


1  -  ( 1-P2 )h2 vi v2 


i  ,  ,  2,2,422  ,22  ,22  „  2 ,  2 

2  (  1  - P  )  h  v  v  +  h  v  +  h  v  +  2P  h  v  i  v  +  1 


1  +  (1-P2)h2viv2 

,  ,  2,2,422  ,22  ,22  „  2  ,  2 

(1-P  )  h  v  v  +  h  v  +  h  v  -  2P  h  v  v  4  1 
11  1  2  12 


(B-22I 


!  V 


k** 


i 

a 

I 

nV 

$ 

sS 

V 


•TjJ 


I 


& 

$ 

Wl 


$ 


'Si 


» 


>:■ 

*<• 

% 


W 


«-*  A.Tjw'ViU.' 


When  v  = 


' 2~  X 


Eq  (B-22)  becomes 


1  +  ( 1-p2 )h2v“ 


,  2,2,44.  2  2  .  2,  2  2 

[_  (  1-P  )  h  v  +  2h  V  -  2P  h  v  +  1 


2  2  2 

1  -  ( 1-P  )h  V 


,  ,  2  ,  2,  4  4  2  2  .  2,  2  2 

(l-P)hv  +  2hv  +  2  P  h  v  +  1 


IB-23) 


Substituting  in  the  extreme  limits  of  P  one  finds  that 


1 

P=0 

,  .2  2 

1  +  h  v 

( B— 2  4  ) 


P=1 


1  +  2  h  2  v  2 
1  +  4  h  2  v  2 


(B-25) 


which  aids  in  confirming  the  validity  for  the  following  reasons.  First, 
at  P=0  the  exponential  jpdf  reduces  to  the  square  of  marginal 
exponential  pdf  (see  Eq  (89)).  Therefore,  the  characteristic  function  at 
p=0  becomes  the  square  of  the  the  marginal  exponential  characteristic 
function  (Eq  99).  As  P  gets  larger  the  second  part  of  the  joint  pif 
becomes  increasingly  significant.  Therefore,  this  first  series,  which 
represents  the  fourier  transform  of  only  the  first  part  of  Eq  (2)  should 
only  give  a  port  ion  of  the  total  solution.  Since  the  total  joint 
•'•haraeterist  ie  function  should  give  a  value  of  one  at  P=l,  the  value  of 
only  the  first  part  of  the  total  solution  should  be  something  lees  than 
one.  Equation  (B-25)  obviously  meets  this  condition. 
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Appendix  C 

The  Truncat ion  of  the  Second  Series  of  the 
Joint  Exponential  Characterist ic  Function 

The  purpose  of  this  appendix  is  to  show  the  method  by  which  the 
second  series  of  the  the  characteristic  function  for  the  joint 
exponential  pdf,Eq  ( A—  37),  is  reduced  for  the  purposes  of  numerical 
calculations.  The  series  is  repeated  hear  as  Eq  ( C— 1 ) . 


8p( 1 -P  ) 
2 


C  2 

"n  =  o 


p2n  F2(n  +  k  +  l  +  .5) 

n  +  1 )2  k!  ( 2  n  + 1  +  k ) ! 


s  i  n  2  [  (n  +  i-:+l  +  .  5)^ 


,  ,  ,,,  2,.  ,2  (n+lc  +  l+  .5> 

(1  +  ( ( 1-P  )hv  )  ) 

z 


:c-i ) 


-1 

/*  »» 

.  - 1 

( 1 -p2 ) hv 

Z 

=  tan 

cv 

z 

=  sin 

„  ( 

1 

+  (  (  1-P2  >hv  >2)l/2 

z  ' 

( C-2 ) 


Although  the  ratio  of  the  gamma  function  and  the  factorial  terms  is 
not  equal  to  one  in  this  series,  it  approaches  one  for  values  of  k  n. 

Therefore,  for  k  sufficiently  Large  compared  to  n,  the  ratio  of  the 
gamma  function  to  the  factorial  terms  is  replaced  with  one  for  the 
remaining  *erms.  When  the  ratio  is  assumed  equal  to  one,  the  geometric 
series  Merit  i ty  used  in  Appendix  B  can  be  utilized. 


km 

m 

i 

m 


(iUCUGWAV.' 


w  ■ 


I 


The  determination  of  exactly  where  to  alter  the  series  and  make 


the  approximation  that 

T2 ( n  +  k  + 1  +  . 5  ) 

-  .  1 

k  ’  (  2  n  +  1  +  k  ) 


( C— 3 ) 


is  of  course  arbitrary.  One  basis  upon  which  to  make  such  a  decision  is 
to  decide  how  close  to  one  the  the  left  side  of  Eq  { C— 3 )  should  actually 
be  and  then  determine  the  number  of  terms  required  to  reach  that  point . 
Using  a  Stirling’s  formula  for  the  evaluation  of  the  Gamma  and  factorial 
terms  [1:257]  a  program  was  written  that  produced  the  data  in  Table 
(C-I).  Table  C-I  lists  the  value  of  k  required  for  the  ratio  to  reach 
minimum  values  of  .98  and  .99  respectively  for  values  of  n  ranging  from 
0  to  20.  The  values  of  k  were  incremented  in  steps  of  100  and  only  the 
even  numbered  n  values  are  given  for  purposes  of  brevity.  The 
documented  program  is  listed  in  the  software  appendix  under  the  title 
GA.MRAT .  BAS . 


TABLE  C-I 


Number  of  k  Terms  Required  to  Reach  .98  and  .99 


n 

.98  k 

.99  k 

0 

1  00 

1  00 

O 

L. 

4  00 

700 

4 

1000 

2  C00 

e 

2  100 

4000 

8 

3400 

4  7  00 

10 

4600 

6  9  00 

12 

6  100 

8900 

14 

7700 

1  1300 

16 

9  100 

12  3  0  0 

18 

1  1200 

1  3  800 

20 

1  3900 

15900 

The  series  can  now  be  summed  up  to  as  many  k  terms  as  required  and 
then  a  remainder  term,  calculated  by  setting  the  left  side  of  Eq  (C-3) 
equal  to  one,  should  be  added.  This  remainder  term  is  found  as  follows. 
After  summing  m  terms,  and  using  Eq  (C-3),  the  series  on  k  is 


p2k  s in  2 [ { n  +  k+ 1  +  . 5 ) 9  ] 

,  , _  2  ,  2  ( n+k+ 1 t . 5) 

(  1  +  (  (  1  -  P  )hv)  ) 


(C-4 ) 


To  use  the  geometric  expansion  identity,  the  lower  index  of  the 
summation  must  be  zero.  This  is  accomplished  by  letting  q  =  k-m  and 
thus  k  =  q+m.  Then  Eq  (C-4)  becomes 


4 m  e  i  2q  ,  2  r  .  „  _  _  , 

P  \  P  sm  [  ( n  +  q  +  m+ 1  +  .  o  )  0  ] 

2  Cn  +  m+1.5)/  2  q 

(  l  +  r  )  (i+/)4 


(  l  +  r2  )q 


2  2 

where  i  -  ((1-P  )hv)  .  Dealing  only  with  the  summation  now,  Euler’s 
identity  can  be  used  to  expand  the  sine  term,  the  series  thus  becomes 


K-> 


p2q(exp(j(n+q+m+l+.5)£)-exp(-j(n+q+m+l+.5)@))2 

2  2  u 

(  2.7  )  (1+A  ) 


(C-6) 


.*'(1 

y;l 


I*. 


Expanding  Eq  (C-6),  the  summation  can  be  written 


(  1  +  r  ) 


P2e!j9  lq 


(1+/  ) 


exp( 2j0(m+n+l . 5 ) ) 


2  -  2  j  0  -|  q 
P  e 


(  1  +  /  ) 


exp(-2j0(m+n+1.5) )  ( C— 7 ) 


Letting 


and  a  =  0( m+n+ 1.5) 


<1  +  7  ) 


(C-8 ) 


the  application  og  the  geometric  expansion  identity  brings  the  series 
remainder  to 


-2jO  1 


1-1  1-E 


T  _  2  j  0 


1-Ee 


2  j  0 


( C-9 ) 


From  this  point  on  the  manipulations  are  not  unlike  those  in 
Appendix  B.  Rather  than  presenting  the  algebraic  gobbledygook  the  final 
solution  is  just  presented  in  Eq  (C-10). 


P  1  isin(<7)  Ac  os  (O) 

o,  i  .  2  (n  +  m+  l  5)  ~  T  +  r  2  _  2  ~  ,  2  TT  <  f  ~  1 0  ) 

2  (  1  +  I  )  1  +i  A  +  i  A  +  £ 


i  +  ; 2 ) 


a  =  0( m+n+ 1.5) 


/  =  ( ( 1-p2 )hv)2 


=  1-E cos (20) 


=  isinl 20) 


This  result  is  used  in  the  program  MAIN2 . BAS  which  computes  values  of 
the  exponential  characteristic  function. 


Appendix  D 


Software 


This  appendix  contains  the  source  code  for  some  of  the  programs 
that  were  written  in  direct  support  of  of  this  thesis.  Those  programs 
written  in  support  of  Part  A  of  the  analysis  chapter  are  written 
primarily  in  BASIC  and  were  run  under  compiled  TurboBasic  .  Those 
programs  written  in  support  of  Part  B  are  in  Fortran  and  were  compiled 
and  run  on  the  Vax  11/785  that  was  provided  by  AFW.AL.  Each  program  has 
its  own  documentation  and  abstract;  the  titles  of  the  programs  are 
listed  below. 

EKPOPDF.BAS:  Calculates  values  for  the  exponential  .joint 
probability  density  function. 

MAIN. BAS:  Computes  values  for  the  .joint  exponential 
characteristic  function. 

GAMRAT . BAS :  Determines  the  ratio  of  the  gamma  and  factorial 
terms  in  the  second  series  of  the  exponential  joint  CF. 

RC32.BAS:  Performs  the  numerical  integration  to  get  the 
<RCS>/Area  for  exponential  surfaces. 

CORCOEFF . FOR :  Builds  a  correlation  matrix. 

CHOLESKY . FOR :  Finds  the  Cholesky  square  root  of  the 
correlation  matrix. 

GES1  "RF . FOR :  Generates  random  rough  surfaces  for  different 
statistics. 

Bl.TLD.FOR:  Buils  the  data  file  for  use  by  the  RCS-Basio 


Scattering  Code  program. 


10  sud  EXPOPDF( v, p, zl , z2 , w)  static 
•  *********************************************************** 

’  *  * 

’*  Sub:  EXPOPDF . BAS  * 

’  *  * 

’*  Purpose:  Calculates  the  value  of  the  exponential  * 

’*  pdf  proposed  by  Dr.  Pyati.  * 

’*  * 

’ *  Input  Parameters :  * 

’*  h  =  RMS  surface  height  * 

’*  p  =  correlation  coefficient  * 

’*  zl  =  height  at  a  point  1  * 

’*  z2  =  height  at  a  point  2  * 

’  *  * 

’ *  Output  Parameters :  * 

'*  w  =  probabilty  of  getting  z2  and  zl  * 

’  *  * 

>  *^***:M** ************************************************* 
pi  =  4*atn( 1 ) 

’  Calculate  values  for  the  sugnum  functions  and  take  absolute 
’  values  of  zl  and  z2 
sgnzl  =  1  :  sgnz2  =  1 
if  zl  <  0  then  sgnzl  =  -1  end  if 

if  zl  =  0  then  sgnzl  =  0  end  if 

if  z2  <  0  then  sgnz2  =  -1  end  if 

if  z2  =  0  then  sgnz2  =  0  end  if 

zl  =  abs(zl)  :  z2  =  abs(z2) 

’  Calculate  the  first  term 

terml  =  e\p( -( zl+z2 ) / (h* ( l-p~2 ) ) )/( 4*h*2* ( l-p~2 ) ) 

’  Compute  the  zeroth  ordered  Bessel  function 
30  a  =  0  :  x  =  sqr(  zl  *z2  )  *2*p/(h*  (  l-p~2  )  )  :  xx  =  \' 

40  call  modbesl (a,xx , term2 )  ’  the  modified  bessel  function 

’  Sum  only  the  first  15  terms  of  the  infinite  seires  of  Bessel 
’  functions 
50  sum  =  0 

GO  for  n  =  0  to  15 

70  a  =  2*n+l  :  xx  =  x 

3 0  ca 1 1  modbes 1 ( a , xx , bessy ) 

90  ■  tm  =  bessy/a '2 

100  sum  =  sum+tm 

110  if  sum  :  oldsum  then 

120  goto  500 

1 30  end  i  f 

140  oldsum  =  sum 

150  next  n 

500  '  erm3  =  sum*  (  3/pi~2  )  *sgr.z  1  *sgnz2 

’  fombine  the  terms  to  get  the  total  solution 
510  u-  =  term  1  *  ( term2  +  temi3  ) 
end  sub 


V 
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sub  modbesl(a,.\,bj )  static 

Sub:  MODBESL 

Purpose:  Returns  the  value  of  the  modified  bessel 
function  I 

Parameters : 

a  =  order  of  the  bessel  function 
x  =  point  at  which  bessel  is  to  be  calculated 
bj  =  returned  value  of  the  bessel  function 

i  i  y  i  t  i  i  9  t  t  i  t  i  t  t  i  i  i  f  i  y  J  i  )  y  i  i  i  i  )  r  t  f  f  t  i  y  i  i  *  t  y  i  f  f  i  i  i  i  f  *  i  y  i  t 

if  a  =  0  then  ’  Compute  by  polynomial  approximation 
if  abs(x)  <  3.75  or  abs(x)  =  3.75  then 
t  =  x/3.75 

bj  =  1  +  3.5156229*t~2  +  3 . 0899424*t~4 
bj  =  bj  +  1 . 2067492*t~6  +  . 2659732*t"8 

bj  =  bj  +  .0360768*t~10  +  .0045813*1/12 

goto  done 
end  if 

if  x  >  3.75  then 
t  =  x/3.75 

bj  =  .3894228+  .01328592/t  +  .00225319/t~2  - 
bj  =  bj  - . 00157565/t^3  +  .00916281/t~4 
bj  =  bj  -  .02057706/^5  +  .02635537/t~6 

bj  =  bj  -  .01647633/t~7  +  . 00392377/t~8 

bj  =  bj/ ( sqr(x ) *exp{ -x ) ) 
goto  done 
end  if 
end  if 

if  a  =  1  then  ’  Compute  by  polynomial  approximation 
if  abs(x)  <  3.75  or  abs(x)  =  3.75  then 
t  =  x/3.75 

bj  =  1/2  +  .87890594*1/2  +  . 54918869*t~4 
bj  =  bj  +.15084934*tA6  +  . 02658733*t~8 
bj  =  bj  +  .00301532*1/10  +  .  0003241 1  *t~  12 
bj  =  x*bj 
goto  done 
end  if 

if  x  >  3.75  then 
t  =  x/3.75 

bj  =  .39894228  -  .3988024/t  -  .00362018/1/2 
bj  =  bj  +  . 00163801/t~3  -  .01031555/t~4 

bj  -  bj  +  .02282967/1/5  -  . 02895312/t~6 

bj  -  bj  +  .01787654/t~7  -  . 00420059/t~8 

bj  =  bj/(sqr(x)*exp(-x) ) 
goto  done 
end  i  f 
end  i  f 


’ ( sub  MODBESL . BAS  cont . ) 

’  Use  infinite  series  representation 
for  m  =  0  to  10000 
gosub  NEWTERM 
bj  =  bj  +  term 
if  TERM  <  IE- 15  then 
goto  done 
end  if 
next  m 

NfEWTERM: 

numer  =  (x/2P(2*m  +  a) 
for  n  :  1  to  m 

numer  =  numer  /  n 
next  n 

for  n  =  1  to  (m+a) 

numer  =  numer  /  n 
next  n 

return 

999  done : 

end  sub 
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********************************************************** 


* 

Prgm : GAMRAT . BAS  * 

* 

Purpose:  This  program  calculates  the  ratio  of  the  * 
gamma  function  term  to  the  factorial  terms  in  the  * 
part  of  the  exponential  CF  that  was  not  put  into  * 
closed  form.  .4s  set  up,  this  program  will  output  * 
the  value  of  k  needed  for  the  ratio  to  be  at  least  * 
.99.  This  data  is  then  used  in  the  program  that  * 
calculates  the  CF  (MAIN2.BAS).  * 

* 

Algorithm:  The  program  uses  Stirling’s  formula  * 

[Abramowitz  &  Stegun:  257]  to  represent  the  * 

values  of  the  gamma  functions  and  factorial  * 

terms.  Logaritmic  calculations  are  utilized  to  * 
avoid  overflow  problems.  * 

* 


’ ********************************************************* 
for  n  =  0  to  20  step  1 

for  k  =  0  to  100000  step  100 

ul  =  log(n+k+l . 5 )  :  um  =( l+l/( 12*(n+k+l .5) ) )+< l/(288*(n+k+l .5 P2) ) 
um  -  log ( um) 

vl  =  log( k+1 )  :  v2  =  log ( k+2*n+2 } 

vml  =(l  +  l/(12*(k+l)))  +  l/<288*(k+ir2)  :  vml  =  log(vml) 
vm2  =(  l  +  l/(  12*(2*n+k+2)  )  ) +  1/ ( 288* ( 2*n+k+2 T 2 )  :  vm2  =  log(vm2> 
term  =  ( 2*n+2*k+2 ) *ul+2*um-k* ( vl+v2 )-2*n*v2-vl/2-l . 5*v2-vml-vm2 
if  exp ( term) > .99  then 

print  using  "  n  =  k  =  ******  ratio  - 

**. *********" ;n,k,exp( term) 
goto  100 
end  if 
next  k 
100  next  n 


ttt**ttttt*t**ttt****t**t******tt*t*ttttt*****t*tt*******tttt*t 
t  * 

*  Prgm:  MAIN2.BAS  * 

t  * 

*  Purpose:  Creates  the  values  of  the  joint  characteristic  * 

*  function  of  the  exponential  PDF.  For  correlation  * 

*  coefficient  very  near  1,  this  program  takes  a  long  * 

*  time  to  calculate  a  single  data  point.  As  the  * 

*  correlation  coefficient  decreases,  the  infinte  series  * 

*  in  the  solution  becomes  less  significant  since  it  is  * 

*  multiplied  by  the  correlation  coefficient  raised  to  a  * 

*  high  power.  Therefore,  the  program  speeds  up  as  it  * 

*  progresses.  However,  this  program  still  took  56  hours  * 

*  to  complete  using  compiled  TurboBASIC  and  the  8087  * 

*  numerical  coprocessor  in  the  IBM  PC.  * 

*  * 

*  Algorithm:  .As  with  GAMRAT.BAS,  Stirling’s  formula  is  used  * 

*  for  the  gamma  functions  and  factorial  terms  while  * 

*  logarithmic  calculations  are  used  to  avoid  overflow.  * 

*  Otherwise,  this  program  is  simply  a  straight  forward  * 

*  implementation  of  the  formulas  derived  in  Appendixes  * 

*  A  through  C.  * 

*  * 

*  Input  Files:  h:one  * 

*  Output  Files:  E:\THESIS\EXPO.DAT  * 

*  * 
tttttttttt**t**t*******tt*t***********tttt**t**tttt*tttttttt ttt 

Use  double  precision  variables  to  reduce  error, 
defdbl  a-h,o-z  :  defi  h  i-n 

Setup  parameters  for  the  run  of  5001  data  points  spaced  .0002 
correlation  lengths  apart  with  h*vz  =  10 
pi=4*atn( 1 ) 

v  =  10.0  ’  Value  of  Vz 

space  =  .0002  ’  Spacing  between  data  points 

h=1.0  ’  RMS  value  of  the  surface  heights 


h=1.0  ’  RMS  value  of  the  surface  heights 

npoints  =  5001  ’  Total  number  of  points  to  calculate 

Open  a  random  access  file  to  store  the  data  points  in 
fout $z"e : \thesis\expo.dat” 
open  foutS  as  si  len=8 
field  s] ,  8  as  cvalS 

Increment  from  0  correlation  lengths  to  1  correlation  length  in 
steps  of  .0002  cor.  lens, 
for  ntow  r  0  to  npoints 
tow  =  ntow*space 

p  =  exp(-tow~2)  'Gaussian  Cor  Coefficient 
gosub  char fun  'Calc  the  CF 
'  Provide  feedback  to  the  user 

print  using  "  t-- . ssss  p charts . =  " ; t  ow , p , char 
’  put  the  data  in  the  file 
Iset  ova  1 S=mkdS ( char )  :  put  si 


next  ntow 
close  *1 
end 

500  eharfun: 

yyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyy 

’  Sub  CHARFUN 

y  y 

’  Purpose:  Evaluates  the  sums  of  the  firts  and  second  ’ 

’  series  of  the  characteristic  function  solution  ’ 

y  y 

’  Input  Parameters:  ’ 

'  h  =  RMS  value  of  the  surface  heights  ’ 

’  v  =  Absolute  vale  of  Vz  ’ 

’  p  =  Correlation  Coefficient  ’ 

y  y 

’  Returned  Parameters:  ’ 

'  char  =  Value  of  the  characteristic  function  ’ 

y  y 

yyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyjyyyyyyyyyy 
y 

’  These  are  the  99%  cutoff  values  where  the  ratio 
’  of  the  gamma  and  factorial  terms  becomes  >  .99 
’  ( see  Appendix  C )  and  program  GAMRAT . BAS 

dim  maxk( 20) 

maxk(0)  =100  :  maxkt 1)  =200 

maxk(2)  =400  :  maxk(3)  =600 

maxk(4)  =1000  :  maxk(5)  =1500 

maxk(6)  =2100  :  ma_xk(7)  =2800 

maxk(8)  =3400  :  maxk(9)  =4200 

maxkt 10) =4600  :  maxk( 11 )=5000 

maxk( 12) =6 1 00  :  maxk( 13 ) =6600 

maxkt 14 ) =7700  :  maxkt 15) =8300 

maxkt 16 ) =9100  :  maxk ( 17 )= 10600 

maxkt 18) =11 200:  maxkt 19) =12000 
maxkt  20)  =  13900 

sum2  =  0  ’  Initialize  variable 

if  p  =  0  then  goto  900  ’  Only  the  closed  form  portion  is  relevant 
if  p  =  1  then 
char  =  1 
goto  999 
end  if 

’  The  portion  of  the  solution  which  is  left  in  summation  form 
’  (see  Appendix  C)  is  calculated  first.  The  outer  summation 
’  is  truncated  after  21  terms  whereas  the  inner  series  is 
’  truncated  based  upon  the  ration  of  the  gamma  and  factorial 
’  terms  for  a  given  (n).  A  remainder  term  is  then  added  to 
’  compensate  for  the  truncation  of  the  inner  series. 

’  Set  up  parameters  that  are  constant  within  the  loops 
e  =  ( 1 -p~2 ) *h*v 

t  =  atn(o) 

1>1  =  I<>gt  p)  :  =  logt  !+oA2)  :am=8*  (  1  -p~2  ) /pi  “2 
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for  n  =  0  to  20  ’  Truncate  the  series  at  20 

anl  =  2*log( 2*n+l ) 
termk  -  0 

for  k  =  0  to  (maxk(n)-l) 
ul  =  log(n+k+1.5) 
um  =  log( 1+1/ ( 12 . 0* ( n+k+1 . 5 ) ) ) 
vl  r  log(k+l) 
v2  =  log(k+2+2 .0*n) 
vml  =  log(l+l/(12.0*(k+l) ) ) 
vm2  =  log( 1+1/ ( 12.0*<2.0*n+k+2) ) ) 

U— ( 2 • 0*n+2 . 0*k+2 ) *ul+2*uml - ( k+ . 5 ) *vl 

u=u- (k+2 . 0*n+l . 5 ) *v2-vml-vm2 

st  =  (sin( (n+k+1 . 5 ) *t) ) ~2 

if  st  =  0  then  goto  790 

tnum  =  pi* ( 2 . 0*n+2 . 0*k+l ) +u 

tden  =  anl+cl*  ( n+k+1 . 5 ) 

term  =  am*exp( tnum- tden ) 

if  term  <  IE-40  then  goto  795  ’  to  small  to  matter 
term  =  term*st 
termk  -  termk  +  term 
next  k 

gosub  remainder  ’  Add  in  the  compensation  for  the  truncation 
sum2  =  sum2  +  termk  +  remain 
next  n 

gosub  sumone  ’  Get  the  closed  form  part  of  the  solution 

char  =  suml  +  sum2 

return 

sumone : 
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’  Sub:  SUMONE 
} 

’  Purpose:  Calculate  the  first  series  which  has 

’  a  closed  form  solution  (see  Appendix  B). 

» 

'  Algorithm:  Implements  Eq  B.24 
» 

’  Input  Parameters: 

h  =  RMS  value  of  the  surface  heights 
’  v  =  Vz 

'  p  =  Correlation  Coefficient 

» 

'  Output  Paramteres: 

'  suml  =  Summation  of  the  first  series 

si  =  1  +  h~2*v~2*( l+p~2 ) 

s!2  =  l-h~6*v~6* ( p~6-3*p" !+3*p~2-l ) 

si  2  =  sl2-h~4*V'4*(p~4+2*p'2-3)+l-r2*v~2*(p~2+3) 

suml=  si /si  2 

return 


fmm 


1200  remainder: 

>  >  »  t  »»>»»»>»  f  I  I  t  71  t  t  t  ft  t  »»»>♦>  t  f  1  J  t  t  t  t  f  1  f  t  t  ft  1  t  1 

i  t 

’  Sub:  REMAINDER 

)  f 

’  Purpose:  The  second  series  of  the  solution  has  ’ 

’  been  truncated  at  the  point  wheie  the  ratio  of  ’ 

’  the  Gamma  functions  and  the  factorial  terms  is  ’ 

’  at  least  .99.  This  routine  calculates  the  ’ 

’  remainder  of  that  term  of  the  series  with  the  ’ 

’  assumption  that  the  gamma  and  factorial  terms  ' 

’  completely  concel .  ’ 

i  y 

1  Algorithm:  Implements  Eq  C.10  ’ 

y  y 

’  Input  Parameters:  ’ 

’  h  =  RMS  value  of  the  surface  heights  ’ 

’  v  =  Vz  ’ 

’  p  =  Correlation  Coefficient  ’ 

’  t  =  theta  angle  ( see  Sub  CharFun )  ' 

’  c  =  shorthand  ( see  Sub  CharFun )  ’ 

’  maxk(n)  and  n  =  (see  Sub  CharFun)  ’ 

’  pi,  anl  and  cl=  (see  Sub  CharFun)  ’ 

i  y 

’  Output  Paramteres:  ’ 

’  remain  =  the  remainder  ’ 

f  » 

yyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyy 

y 

sigma  =  2* (maxk(n)+n+l . 5 ) *t 

a  =  l-p~2*cos(2*t)/( l+c"2)  :  al  =  log(a) 

b  =  p~2*sin( 2*t ) / ( l+e~2 )  :  bl  =  log(b) 

albl  =  log(a"'2+b"2 ) 

el  =  log( l-p'2/( l+c*2) ) 

pl2  =  log( l-p"2 ) 

rl  =  pi*  (  2*  (maxk( n ) +n) +  1 ) -anl+pl2-cl*  ( maxli ( n ) +n+l .  5  ) 

rml  =  rl-el 

rm2  =  rl+bl-albl 

rm3  =  rl+al-albl 

remain^  ( exp(  rml )  +sin(  sigma  )  *exp(  rm2  )cos(  sigma )  *e\-p(  m3  )  ) 

remain  -  remain*4/pi '2 

return 


************************************************************** 
*  * 

*  Prgm:  RCS2.BAS  * 

*  * 

*  Purpose:  Calculates  the  values  of  the  normalized  * 

*  RCS  of  the  exponentially  rough  surface.  * 

*  * 

*  .Algorithm:  Integrates  Eq  (87)  using  a  Simpson’s  * 

1  *  approximation  routine  and  the  5001  data  points  of  the  * 

'*  eponential  CF  generated  by  MAIN2.BAS  * 

'*  The  Bessel  function  is  calculated  using  a  polynomial  * 

'*  approximation  taken  from  [Abramowitz  &  Stegun:369}  * 

1  *  * 

'*  Input  Files:  E:\THESIS\EXPO.DAT  * 

1  *  * 

’*  Output  Files:  E:\THESIS\JEXPO.DAT  * 

’  *  * 

'************************************************************** 


defdbl  a-h,o-z 
defint  i-n 

'  expo  will  hold  the  data  points  for  the  exponential  CF 
’  value  will  hold  the  integrand  value  at  the  5001  discrete 
’  points 

’  result  will  hold  the  vaule  of  the  integral  for  the  15 
’  different  angles  and  the  6  values  of  atan(s) 
dim  value) 5001 ) , result) 15,6) ,expo( 5001 ) 
pi=4*atn( 1 ) 

’  This  section  takes  care  of  the  initialization  of  the 
'  arrays  that  are  used  in  the  bessel  function 
’  subroutine  rather  than  reintial izing  the  arrays 
'  each  of  the  5001  times  the  routine  is  called, 
dim  ja# ( 7 ) , jb# ( 7 ) , jo* ( 7 ) 

’  Data  for  calcualtion  of  Bessel(x)  for  x<=3 


ja# ( 0 ) 
ja# ( 1 ) 
ja=( 3 ) 


=  1.00 

=  -2.2499997 
=  -.3163866 


ja#(5 )  =  -.0039444 


ja* ( 2 )  =  1 . 2656208 
ja*<4)  =  .044447900 
ja*(S)  =  .0002100 


Data  for  calcualtion  of  Bessel(x)  for  x>3 


jo=(0)  = 
jc=(2)  =  ■ 
jo* ( 4 )  = 
jc== ( 6 )  = 


. 79788456 
■ . 00552740 
.00137237 
.00014476 


jc#( 1 ) 
jc*(  3 ) 


-.00000077 

-.00009512 


jc*(5)  =  -.00072805 


jb# ( 0 )  =  3.00 
.jb* (  2  )  =  -.04166397 
jb# ( 4 )  =  .00262573 
jb# ( 6 )  =  -.00029333 


jb# ! 1 )  = 


'8539816 


jb# ( 3 )  =  -.00003954 
jb# ( 5 )  =  -.00054125 
jb# ( 7 )  =  .00013558 


open  "E:\THESIS\EXPO.DAT"  as  *1  len=8 
fi°ld  *1,  8  as  cvalS 
print  "Getting  data  from  data  file." 
fur  n=l  to  5001 


tk  nn  wTtwv  fvwv  wv  rjrjfvn'  rv  fv 


get  #1  :  expo(n)=cvd(cval$) 
next  n 
close  #1 
(RCS2.BAS  cont. ) 

spacer . 0002 
h=l .0: vz=10.0 

xl  =  1/ (  l+h^2*vz,'2 )  ’  One  dimensional  Char  Func 

for  ang%  =  0%  to  70%  step  5% 

nang  =  nang+1%  ’  counter  for  result (nang, ns lope) 
rang  =  pi*ang%/180  ’ 

for  m  =  5%  to  30%  step  5% 

nslope=nslope+l  ’  Slope  index  counter 
’  Scale  the  wavenumber  (ck)  according  to  the  angle  of 
’  incidence  and  teh  mean  surface  slope 
ck  =  10/ ( tan(m*pi/180 ) *cos ( rang ) ) 
vx  =  2*ck*sin( rang !  ’ 

for  n  =  1  to  5001 
towr ( n-1 ) * space 

arg  =  vx*tow  *  Bessel  argument 
if  (arg<3  or  arg=3)  then 
gosub  BESLLOW 
else 

gosub  BESLHI 
end  if 

vaiue(n)  =  tow*  j0#»  ( ex-po(n)-xl  '2 )  *2*ck~2 
next  n 

locate  1,1 ‘.print  ang%,m 
npoints=5001 

call  simpson ( %alue( ) .npoints .space, result! nang, nslope) ) 
next  m 
next  ang% 

’  Transfer  data  to  a  file 

open  "E:\THESIS\JEXPO.DAT"  for  output  as  ?2 
for  nslope=l  to  6 

for  nang  =  1  to  15 

print  42 ,nslope*5 . 0 , (nangl ) *5 . 0 , result ( nang , nslope ) 
next  nang 
next  nslope 
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!  Zeroth  order  Bessel  function  polynomial 

approximations.  Taken  from  Abromowitz  and  Stegun. 

’  Input  Parameters : 

'  Coefficient  Arrays  =  Filled  in  main  prgm 
arg  =  argument  of  the  bessel  function 

’  Returned  Parameters: 

’  jO#  =  Value  of  the  bessel  function 

^tiyitifitiiitfiiyflyyyyyyyyyyyyyyyyyyyyyyyyyyjyyjyfffy 

BESLLOW: 

j0#=0.0 

for  nbl%  =  0  to  6 

jO#  =  jO#  +  ja#(nbl%)*(arg/3)A (2*nbl%) 
next  nbl% 
return 

BESLHI: 

theta#=0.0:f0#=0.0 
for  nbl%  =  0  to  7 

theta#  =  theta#  +  jb#( nbl%) * ( 3/arg ) * ( nbl%-l ) 
next  nbl% 

for  nbl%  =  0  to  6 

fO#  =  fO#  +  jc#( nbl%) * ( 3/arg ) "nbl% 
next  nbl% 

jO#  =  argA ( -1/2 ) *fO#*cos( theta# ) 
return 

sub  simpson(f#( 1) ,n%,h#, result#)  static 

tiyiyiiyiitiyyyyyyyyyyyyytyy  yyyyyyyyyyyyyyyyyyyyyyyyyyy 


Sub:  SIMPSON' 

Purpose:  Performs  a  Simpson's  approxomation 
Integration.  Algorithm  taken  from  'Applied 
Numerical  .Analysis  by  Curtis  F.  Gerald. 

Input.  Parameters: 

f#! 1 )  =  a  one  dimensional  array  of  data  points 
n%  s  the  number  of  data  points  in  the  array 
hs  =  the  spacing  between  the  data  points 

Returned  Parameters: 

results  =  the  approximated  value  of  the  integral 
local  i%, npanel%, nbegin%, nhal f% 


Check  to  see  if  the  number  of  panels  is  even 
npanel%  =  n%  -  1 
nhalf%  =  int(npanel%/2 ) 
nbegin%  =  1 
result%  =  0 

if  npanel%  -  2*nhalf%  =  0  then 
goto  EVEN 
end  if 

’  Number  of  panels  is  odd.  Use  3/8  rule  on  first  three 
’  panels  and  1/3  rule  on  the  rest  of  the  panels, 
result*  =  3*h#/8*  (  f#<l)  +  3*f#(2)  +  3*f#(3)  +  f*(4)  ) 
nbegin%  =  4 


EN :  ’ 

’  Apply  1/3  rule,  add  in  first  , second,  and  last 
’  values. 

rtemp#=  h#/3  *  (  f#(nbegin)  +  4*f#(nbegin+l )  +  f#(n)  ) 
result#  =  result#  +  rtemp# 
nbegin%  =  nbegin#  +  2 
if  nbegin%  =  n%  then 
goto  DONE 
end  if 

’  The  pattern  after  nbegin+2  is  repetitive. 

( RCS2 . BAS  subprogram  SIMPSON  cont . ) 

’  Get  nend,  the  place  to  stop. 
nend%  =  n%-2 

for  i%  =  nbegin%  to  nend%  step  2 

result#  =  result#  +  h#/3  *  (  2*f#(i)  +  4*f#(i+l)  ) 
next  i% 


IE: 

end  sub 


C* ********************************************************** 


c  * 

C  Prgm:  CorCoeff.FOR  * 
C  * 
C  Purpose:  Fortran  program  to  create  a  file  containing  * 
C  the  values  of  the  correlation  matrix.  * 
C  * 
C  Algoritm:  Corrleation  is  a  function  of  the  distance  * 
C  (d)  between  two  points.  (d)  is  calculated  by  finding  * 
C  the  horizontal  (b)  and  the  vertical  (a)  distance  * 
C  between  two  points  and  then  solving  the  right  * 
C  triangel  equation  to  find  the  hypotenues  (d).  * 
C  Following  this,  the  Gaussian  Correlation  coefficient  * 
C  equation  is  applied.  * 
C  * 


C* ********************************************************** 

parameter) n=15 , n2=225 , space= . 3 ) 
virtual  coeff(n2,n2) 
real* 16  coeff,a,b,d 
integer  corrow, cored , ref row, refcol 

C  Create  lower  triangle  of  corr.  matrix  (symetric  matrix) 

do  33  i=l ,n2 

C  calculate  the  row  point  (i)  is  in  (reference  Fig. (21)) 

corrow  =  ( { ( i-1 )/n)+l ) 

C  claculate  the  column  point  (i)  is  in 

cored  =  i-(  ( corrow- 1 )  *n) 
write(6,*)  i 
do  31  j=l,i 

C  calculate  the  row  point  (.;)  is  in 

refrow  =  ( ( ( j-1  )/n )  +  l ) 

C  calculate  the  column  point  i.j)  is  in 

refcol  =  j  -  (( refrow- 1 ) *n) 

C  the  vertical  distance  (a) 

a  =  ( corrow- ref row )* space 

C  the  horizontal  distance  (b) 

b  =  ( corcol-refcol ) *space 

d  =  sqrt ( a**2+b* *2 ) 

eoeff(i,j)  =  exp(-d**2) 

31  continue 

33  enrit  inue 

■  >[^en(  uni  t=  1 1 ,  name:  ’corg3.dat’  ,  status^  ’NEK’  ) 

Ho  43  i  =  l,n2 
write(6,*)  i 
do  41  j  = 1 ,  i 

write(ll,*)  coeff(i,,j) 

4 1  con  t i nue 

43  continue 

close! unit=l 1 ) 


C************************************************************** 
C  * 

C  Prgm : CHOLESKY . FOR  * 

C  * 

C  Purpose:  Finds  the  Cholesky  square  root  of  the  * 

C  correlation  matrix.  * 

C  * 

C  Algorithm:  Implements  Eqs  (xx>  through  ( xx ) .  * 

C  * 

C  Input  Files:  The  correlation  matrix  file  created  by  * 

C  CORCOEFF . FOR  * 

C  Output  Files:  The  lower  half  of  the  Cholesky  matrix.  * 

C  * 

C************?************************************************* 

parameter! n2=225 ) 

virtual  cor ( n2 , n2 ) , cho ( n2 , n2 ) 

real *  16  cor, cho, sum 

open(unit=ll ,name=’ corg3.dat’ ,status=’OLD’ ) 

C  Input  the  correlation  matrix  (lower  triangle) 

do  33  i=l,n2 
wt  i te ( 6 ,  *  )  i 
do  31  j= 1 , i 

read (11,*)  cor(i,j) 

31  continue 

33  continue 

close! unit=l 1 ) 

C  decompose  into  Cholesky  sqrt . 

do  43  k=l ,n2 

do  41  i=l , (k-1 ) 
sum  =  0.0 
do  40  j  =  1 , ( i  —  1 ) 

sum  =  sum+cho(  i  ,  ,j  )  *cho(  k,  j  ) 

40  continue 

cho ( k , i ) = ( cor ( k , i ) -sum ) /cho ( i , i ) 

41  continue 
sum  =  0.0 

do  42  j=l , (k-1 ) 

sum  =  sum+cho  (  k  ,  .j  )  **2 

42  continue 

cho'k.k)  =  sqrt ( cor! k , k ) -sum) 

43  continue 

C  Ouput  data  to  a  file 

open( uni t- 1 1 , name= ’chog3.dat ’ , status = ’ NEW* ) 
do  53  k= 1 , n2 

do  5 1  i  =  1 , k 

write ( 11,*)  cho! k , i  ) 

5 1  con  t i nue 

wri te ( 6  ,  *  )  k 
53  continue 

close ( un  L  t  —  1 1 ) 

stop 

end 


I 
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c* ************************************************************* 
c  * 

C  Prgm:  GESURF . FOR  * 

C  * 

C  Purpose:  Generate  the  data  which  represents  the  rough  * 

C  surface.  * 

C  * 

C  Algorithm:  Implements  Eq  (xx)  to  generate  a  set  of  * 

C  correlated  data  points,  each  of  which  represents  the  * 

C  surface  height  at  a  given  point.  There  are  two  basic  * 

C  parts  of  this  program.  The  first  part  generates  the  * 

C  vectors  of  uncorrelated  Gaussian  and  exponential  surface  * 

C  heights  using  a  Monte-Carlo  technique.  The  second  part  * 

C  multiplies  this  vector  by  the  Cholesky  sqrt  of  the  * 

C  correlation  matrix  to  get  the  vector  of  correlated  * 

C  surface  heights.  * 

C  * 

C  Input  Files:  The  Cholesky  sqrt  file  generated  by  CHOLESKY  * 

C  Output  Files:  A  Gaussian  and  an  exponential  surface  file.  * 
C  * 

C* ************************************************************* 

parameter ( n2-225 ) 

virtual  cho(n2,n2) 

dimension  g ( n2 ) , e ( n2 ) . sg ( n2 ) , se ( n2 ) 

integer*-!  seed 

real *8  ul ,u2 ,api ,gval ,eval 

character*ll  infile, gfile.efile 

wri te (6,10) 

10  format! lx, ’ inputting  cholesky  data  file  ’) 

11  format! lx, ’ input  output  gaussian  file  name’) 

12  format! lx, ’ input  output  exponential  file  name’) 

C  get  data  from  the  Cholesky  sqrt  file. 

open! unit=l 1 , name= ’chog3 .dat ’ ,status= ’OLD ’ ) 
do  34  i=l,n2 


do  33  j  =  1 , i 
read (11,*' 
continue 


cho ( i , j ) 


continue 
close(unit=ll ) 
continue 

Assign  name  for  the  Gaussian  surface  file 
write! 6,11) 
read (5,*)  gfile 

Assign  name  for  the  exponential  surface  file 
write! 6 , 12 ) 
read ( 5 , *  )  e  f  i  1  e 
write! 6 , 40 ) 

Input  a  random  number  generator  seed  to  make  this  a  unique 

set  of  surfaces 

format! lx, ’ input  seed’ ) 

read ( 5 , *  )  seed 

am  =  l*a  tan!  1.0) 


,*Ti’ 


( Prgm  GESURF . FOR  cont . ) 


Monte-Carlo  methos  of  generating  uncorrelated  Gaussian 
Z  surface  heights.  The  maximum  height  possible  of  being 

?  generated  is  5  standard  deviations.  The  variance  is  set 

:  to  1 

50  ul  =  ran(seed) 

ul  =  (ul-.5)*10.0 
u2  =  ran ( seed ) 

gval=exp( -ul**2/2 .0)/sqrt( 2.0*api  ) 
if  (u2  .It.  gval )  then 
g(kount)  =  ul 
kount  =  kount+1 
end  if 

if  (kount  .le.  n2)  go  to  50 

2  Create  the  uncorrelated  exponential  data  by  Monte-Carlo 

kount  =  1 
ul  =  ran  (seed) 
ul  =  ( ul- . 5 ) *  10 . 0 
u2  =  ran (seed) 

eval  =  exp( -sqrt(ul**2 ) )/2.0 
if  (u2  .It.  eval)  then 
e( kount)  =  ul 
kount  =  kount+1 
end  if 

if  (kount  .le.  n2 )  go  to  60 
open  the  output  files 
open(unit=ll .name^gfile.statusr’MEW’ ) 
open ( unit=  12 , name=ef  ile ,  status^  ’  NTIW ’  ) 

perform  the  matric  /vector  multiplication  and  output  data 
do  74  i  =  1 , n2 
sg( i )  =  0.0 
se( i )  =  0.0 
do  73  j  =  1 ,  i 

sg(i)  =  sg( i )+g( j ) «cho( i , j  ) 
se(i)  =  se( i )+e( j ) *cho( i , j ) 
write!  6,*)  cho(i,.j) 
continue 

wTite( 11 , * )  sg(  i ) 
write) 12,*)  se(  i ) 
cont inue 
close(unit=l  1 ) 
close  ( \ini  t=12 ) 

lciop  back  to  create  another  surface  if  desired, 
format ( '  Input  a  1  if  go  again,  -1  to  stop.’) 
wr i te ( 6 , 80 ) 
read (6,*)  nagain 
if  (nagain  .EQ.  1)  then 
go  to  35 
end  if 


MM  >■ 


C* ************************************************************* 

c  * 

C  Prgm:  BUILD. FOR  * 

C  * 

C  Purpose:  Builds  a  data  file  on  the  format  required  by  * 

C  the  Radar  Cross  Section-Basic  Scattering  Code  from  * 

C  the  surface  files  created  by  GESURF.  * 

C  * 

C  Algorithm:  Outputs  the  proper  control  codes  and  * 

C  parameters  to  have  the  RCS-BSC  operate  in  the  perfect  * 

C  conductor,  physical  optics  mode  and  at  the  proper  * 

C  frequency.  .Also,  outputs  the  scaled  surface  data  as  * 

C  392  flat  triangle  so  that  there  are  no  convex/concave  * 

C  surfaces .  * 

C  * 

C  Input  files:  A  surface  file  created  by  GESURF  * 

C  Output  files:  A  data  file  for  input  to  RCS-BSC  * 

C  * 

C************************************************************** 

parameter ( nblockr 1 4 ) 
dimension  z( 15, 15) ,plate( 392 , 3 , 3 ) 
spacer  1 . 0/ ( nblock* 1.0) 
api=4*atan( 1.0) 

C  Input  the  mean  surface  slope  required, 

write (6 , 4 ) 
read  15,*)  s 

4  format ( ’  Input  the  atn(s)  angle’) 

C  Calculate  the  scale  factor  to  multiply  the  heights  by. 

C  Use  the  equation  atan(s)=2h/l  and  the  fact  that  the 

C  surfaces  were  generated  with  h=l.  Thus  h=atan( s ) *1/ ( 2 ) . 

C  Since  there  are  a  total  of  (14*. 3)  correlation  lenghts  in 

C  each  direction  of  the  surface,  a  single  correlation  length 

C  is  equal  to  1/(14*. 3)  units.  Therefore,  the  scale  factor 

C  becomes 

scale= tan ( s*api/ 180 .0) / ( 28 . * . 3 ) 

C  Open  a  surface  data  file  and  read  in  data. 

5  format ( ’  Opened  .dat’ ) 

open ( unitr 1 1 ,namer 'G8.dat' , statusr’OLD’ ) 
wr i te ( 6 , 5 ) 
do  30  irl,l5 
do  20  j  =  1 , 1 5 

read(ll,*)  z(i,j) 
c  write (6,*)  z(i,j) 

20  continue 

30  continue 

close! uni tr 1 1 ) 


o  o 


( Prgm  BUILD. FOR  cont.) 

Calculate  the  x,y  and  z  coordinates  of  each  corner  of  the 
triangle  to  be  input  to  the  program  as  an  individual  plate, 
do  50  n=l,nblock 
write(6,*)  n 
do  40  m=l,nblock 

nop=(n-l ) *2*nblock+ ( 2*m)-l 
plate ( nop , 1 , 1 ) = ( m- 1 ) *space 
plate (nop, 1 , 2 ) = ( n— 1 ) *space 
plate ( nop, 1 , 3 ) =z ( n ,m) * scale 
plate( nop, 2 , 1 ) =m*space 
plate(nop,2,2)=(n-l )*space 
plate ( nop, 2 , 3 ) =z ( n ,m+l ) * scale 
plate ( nop ,3,1) =m*space 
plate ( nop, 3 , 2 ) =n*space 
plate ( nop, 3 , 3 ) =z ( n+1 ,m+ 1 ) * scale 

C  Calculate  total  surface  area  of  the  surface 

delx=(plate(nop,2, 1 ) -plate (nop, 1,1) ) *  *  2 
del z=( plate (nop, 2 , 3 ) -plate (nop, 1,3) ) **2 
xlen=sqrt ( delx+delz ) 

dely= ( plate ( nop ,2,2) -plat e ( nop , 3 , 2 ) ) *  *  2 
delz=( plate (nop, 2 , 3 ) -plate { nop, 3 , 3 ) ) **2 
ylen=sqrt ( dely+delz ) 
area=area+ ( xlen*ylen ) /2 
nop-nop+ 2 

plate(nop,2, 1 )=(m-l )*space 
plate ( nop , 2 , 2 ) =n* space 
plate (nop, 2 , 3 )=z(n+l ,m) * scale 
do  39  nk=l,3 

plate ( nop , 1 , nk ) =plate ( ( nop- 1 ) , 1 , nk ) 
plate ( nop , 3 , nk ) =plate { ( nop- 1 ) , 3 , nk ) 

39  continue 

delx=( plate ( nop, 2 , 1 ) -plate (nop, 3 , 1 ) )**2 
delz= ( plate ( nop ,2,3) -plate ( nop , 2 , 3 ) ) *  *2 
xlen=sqrt (delx+delz) 

dely=( plate (nop, 1 ,2)-plate(nop,2,2) ) **2 
delz=( plate (nop, 1 , 3 ) -plate (nop, 2 , 3 ) )**2 
ylen=sqrt (dely+delz ) 
area=area+ ( xlen*ylen ) /2 

40  continue 

50  continue 

write (6, 75)  area 

75  format  (  ’  .Area  =  *  ,  f  1 0 . 7  ) 

C  Output  to  RCS-BSC  file 

C  while  setting  parameters 

open(  unit-1 1  ,name=  ’G8RCS .  DAT’  , status=  ’NEW ’  ,carr i ageeont- rol  =  ’ '.'CVF ' 
100  format (  f 10.8, ’ , ’ , fl0.8, ’ , ’ fl0.8) 

105  format ( ’ PG: ’ ) 

107  format! ’ 3 ’ ) 

1  1 0  f  o  rma  t  (  ’  TO :  ’  ) 
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115 

format! ’F,F,F,F’ ) 

120 

format! ’F,F,F’ ) 

125 

format! *  F  * ) 

( Prgm  BUILD. FOR  cont.) 

126 

format! ’T,F,F,F,F,F’ ) 

127 

format! ’F,F,F’ ) 

128 

format! ’F,F,F,F’ ) 

C 

(Prgm  BUILD. FOR  cont. 

130 

format! ’UN: ’ ) 

135 

format ( ’ 1 ' ) 

140 

format! ’PD: ’ ) 

145 

format! ’0.0,0.0,90.0,0 

150 

format! ’F,9C.0’ ) 

155 

format! ’-90,90,5’ ) 

160 

format! ’BK: ’ ) 

165 

format! ’0.0’ ) 

170 

format! ’PP: ’ ) 

175 

format! ’T’ ) 

176 

format! ’T,3. ,3. ’ ) 

177 

format! ’ 0 . , 90 . , 10 . ’ ) 

178 

format! ’0. ,30. , 10. ’ ) 

180 

format! ’ LP: ’ ) 

185 

format! ’ T ’ ) 

100 

format! ’ FR: ’ ) 

195 

format! ’5.46’ ) 

2  no 


last  one  is  frequency  Ghz 


write! 1 1 
write! 1 1 
write ( 11 
write ( 11 
write ( 1 1 
write! 11 
write! 1 1 
write! 1 1 
write! 1 1 
write! 11 
write! 1 1 
wr i te ( 11 
write! 1 1 
write! 1 1 
write! 1 1 
write! 1 1 
uri te( 11 
write! 1 1 
write! 1 1 
write! 11 
write! 11 
wr i tell! 
write! 11 
wri  t-"!  1  1 
i  nr i nue 


110) 

115) 

120) 

125) 

126) 

127) 

128) 
130) 
135 ) 
140) 
145) 
150) 
155 ) 
160) 
165) 
170) 

175) 

176) 

177) 
173! 
180) 
185) 
100 ) 
1  05  i 


( Prgm  BUILD . FOR  cont . ) 


Ouput  ( x , y , z )  of  triangle  points 
do  250  i=l , 2*nblock*nblock 
write! 11 , 105) 
write!  11,107) 
do  220  j=l,3 

write! 11 ,100)  plate  (  i  ,.j,  1 )  ,  plate!  i  ,  j  ,  2  )  ,  plate  (  i  ,  j  ,  3 ) 
continue 
continue 
format i ’XQ: ’ ) 
format! ’EN: ' ) 
write! 1 1 , 300 ) 
write! 11 , 305 ) 
close(unit=ll ) 


Appendix  E 


Computer  Generated  Rough  Surfaces  and 


their  Computed  Radar  Cross  Sections 


This  appendix  contains  the  graphs  of  the  computer  generated  Gaussian 
and  exponential  surfaces  along  with  tables  of  the  computed  RCS  values. 
The  tables  contain  the  predicted  RCS’s  for  each  accompanying  figure  at 
discrete  angles  from  normal  incidence  (0°)  to  edge  on  incidence  (90°). 
The  RCS’s  are  tabulated  for  two  different  frequencies,  5.76  Ghz  and 
17  Ghz,  and  for  two  different  values  of  mean  surface  slope  (atan(2h / 1) ) . 
Therefore,  there  are  four  sets  of  RCS’s  values  associated  with  each 
surface.  The  combination  of  a  5.76  Ghz  frequency  and  a  mean  surface 
slope  of  5°  results  in  a  Rayleigh  parameter  of  hk=5.  A  mean  surface 
slope  of  25°  at  5.76  Ghz  results  in  hk=16.  When  the  frequency  is 
increased  to  17  Ghz,  hk=27  at  5°  and  hk=83  at  25°.  Each  table  also 
includes  the  physical  optics  RCS  values  for  a  perfectly  flat  plate  at 
the  two  frequencies  for  comparison.  All  computations  assume  perfect 
conductivity.  Figure  E.21  is  added  to  contrast  surfaces  with  a  mean 
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slopes  of  o  and  2o  . 
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TABLE  E.  1 


<RCS>/Area  in  dBsm 


Gauss ian 


Flat 


aian(s)«5 


16.21 
6.71 
-4 .93 
14.31 
10.70 
-9.72 
-9.15 
12.61 
15.09 
23 . 09 
17.29 
16.43 
17.58 
15.28 
18.57 
31.81 
26.11 
36.35 


hk  »  2T 

19.42 
16.50 
5.84 
11.18 
14.21 
2  5.35 
13.11 
18.38 
15.28 


atanC  s) *25 


11.07 
4.47 
1  .  86 
0 . 59 
7.71 


3 .60 
11.73 


25 . 09 
25 . 06 
16.52 
35.32 
2  1.46 
23.77 
16.06 
13.91 


36 . 20 
10.55 
8.73 
6 .25 
3  .74 
1.18 
-4  .  28 
-11.69 


-99999 


46 . 08 
5 .64 
7.54 
6.13 
-3 . 96 
•20.52 
-5.54 
-7.33 
•20.65 
-5.93 
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■19.59 

■10.08 

■41.25 

■17.77 

■19.83 

•31.13 

■30.60 

■99999 
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TABLE  E . 2 


<RCS>/Area  in  dBsm 


18 

,  34 

16  . 

.  5  5 

19 

.91 

19  . 

.  99 

5 

.  92 

2  . 

,  36 

_  2 

.  6  5 

-  /  , 

.  3  5 

-  5 

.  5  1 

-13  . 

.23 

-3 

.  7  6 

-15. 

.  72 

-18 

.  99 

-  18  . 

.87 

-  1  3 

.03 

-14  , 

,  50 

-  1  7 

.  5  1 

-  18  , 

.  49 

-15 

.  03 

-20  , 

.  6  1 

-12 

.45 

-15. 

.85 

-24 

.  07 

-18  , 

.95 

-15 

.01 

-17  . 

.  9  5 

-23 

.  49 

-3  1. 

,  5  1 

-  16 

.67 

-23  , 

.73 

-19 

.  96 

-24 

.  7  7 

-25 

.  53 

-17. 

.  39 

-25 

.  42 

-33  . 

.  98 

-30 

.  0  3 

-30  . 

.  19 

6  .  35 
0.61 
-4.31 
1  .  86 
6.15 
9 . 80 
7 .68 
5  .  32 
-0.12 
-3.92 
14.56 
17.84 
15.65 
12.46 
2  1.20 
22 . 06 
2  1.03 
15.95 
15.44 


n  <  s  >  *  2  5 

hk  '83 


6 . 29 
4 . 24 


6 . 90 
5 . 46 
3 .64 
-7.68 
5.44 
-0 . 20 
-2.19 
-6.94 
-13.76 
-17.11 
-17.44 
-26.00 
-20.25 
-28 . 97 
-26 .48 
-20 .98 
-15.39 


Flat 

k*  1  1 ♦  k  *  3  5  7 

36  .  20 

46.08 

10.55 

5.64 

8 .73 

7.54 

6 . 25 

6.13 

3.74 

-3 . 96 

1.18 

20.52 

-4  .  28 

-5.54 

-11.69 

-7  .33 

-3.75 

20.65 

1  7  , 

,  30 

-22  . 

,52 

14  . 

.  49 

-19  . 

.  5  9 

-9  , 

,  8  5 

-  10  . 

.  08 

37  , 

.22 

-4  1  . 

.  2  5 

21.47 

28.48 
28.15 
99999 


■19.83 

•31.13 

•30.60 

■99999 


".ABLE  E .  3 


<RCS>/Area  in  dBsm 
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